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Foreword 


This master thesis written by Christian Bodinger is an excellent example of 
the enormous potential multi-sensoral remote sensing has for the derivation 
of biophysical vegetation parameters. It shows that satellite science can help 
us a lot in understanding changes in the environment as a consequence of 
climatic changes. This work, which was done in preparation of the DFG 
funded priority program “EarthShape” to evaluate the impact of vegetation 
on surface formation, received top grades from both reviewers with its scien¬ 
tific approach being close to that of a doctoral thesis. 

The thesis describes a remote-sensing based vegetation classification in four 
different eco zones of Chile, while also investigating the vegetation devel¬ 
opment with high-resolution satellite time-series during El-Nino years and 
drought periods. This is done with approaches from the strongly emerging 
field of machine learning and a critical reflection on the own results, which I 
consider an essential foundation of good scientific work. 

In the meantime, Christian Bodinger has left the Geoinformatics work group 
to join the remote sensing company EOMAP in Seefeld. I am sure he will 
stay open to research questions of modern science and make his way in earth 
observation. I wish him success in doing so and to those interested in the 
topic, a pleasant read. 

Prof. Dr. rer. nat. Volker Hochschild 
(Chair of Geoinformatics, 
Geographical Institute, 
University of Tubingen) 
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Introduction of the Work Group 


The geoinformatics working group at the geographical institute of the Uni¬ 
versity of Tubingen with a 14-year history has its research priorities in the 
field of multi-sensoral remote sensing of hydrological-glaciological process¬ 
es in the high mountainous areas and the recognition of urban structure types 
in high dynamical metropolitan areas. Recently, spectacular results were 
achieved during the monitoring of environmental impacts of refugee camps 
and gravitational mass movement events in subtropical highlands. A regional 
focus is laid on countries of south-east Asia and Africa. Now, the interna¬ 
tional working group has more than 10 employees. 

Besides the remote sensing focus, the working group applies innovative 
methods of analysis (machine learning, etc.) in the field of Geographical 
Information Systems (i.a. Web-GIS-Applications) and uses comprehensive 
geo-relational data bases as well as various numerical modelling approaches 
in hydrology and geomorphology. Worth emphasizing are also the more than 
10-years of cooperation with the University's archeologists for landscape 
reconstruction and the role of culture in early human expansions originating 
from Africa. 

Methods of geoinformatics are practiced during the bachelor’s degree course 
and intensified in the popular master course “Umweltgeographie”. For four 
years the chair also engages in the Tubingen center for scientific advanced 
training with the Geodata Manager qualification, which is also certified by 
the employment agency. 

The working group has comprehensive soft- and hardware suitable for GIS- 
analyses and digital image processing as well as various other equipment 
(field spectrometer, laser particle counter, differential GPS, TDR probes, 
hemispherical camera, UAVs, etc.) to successfully validate remote sensing 
data in the field. Hence, the work group is a popular international coopera¬ 
tion partner for geoscientific projects of all kinds. 

Prof. Dr. rer. nat. Volker Hochschild 
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Abstract 



Graphical abstract to illustrate the workflow behind this thesis. Refer to the Abbrevi¬ 
ations section of this document for clarification of acronyms. 

How is the vegetation distribution influencing the erosion and surface for¬ 
mation in the different eco zones of Chile? To answer this question, it is 
mandatory to possess fundamental knowledge about plant species habitats, 
occurrence and their dynamics. This study utilized both open access Sentinel 
and Landsat satellite imagery as well as the new TanDEM-X Digital Eleva¬ 
tion Model (DEM) in combination with different machine learning environ¬ 
ments (Random Forest, Support Vector Machines). The aims were to derive 
maps of land use and land cover (LULC) in four study sites along a climatic 
gradient and monitor vegetation using monthly Normalized Difference Vege¬ 
tation Index (NDVI) time series. NDVI time series build on catchment-wide 
means in a Mediterranean (La Campana National Park) and a humid- 
temperate (Nahuelbuta National Park) study site were analyzed statistically 
using the BFAST (Breaks For Additive Seasonal and Trend) algorithm 
(Verbesselt et al., 2010) and later on split into class-wise and altitude specific 
series. The mapping outcomes were in support of Random Forest, which 
outperformed the classifications based on Support Vector Machines in each 
of four study sites. Overall Accuracies (OA) for the final maps were consist¬ 
ently above 85% with a coverage of ~10 000 km 2 distinguishing 19 different 
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Abstract 


LULC classes. Time series analysis detected a significant impact of the 
2015/16 El Nino on the NDVI pattern in La Campana but did not find any 
evidence in Nahuelbuta. Vegetation in both sites showed a progressive re¬ 
covery from the 2011 -2015 central Chile megadrought (MD). The findings 
can contribute to a better understanding of climate impacts on Chilean vege¬ 
tation and form a basis of landscape evolution models that include the vege¬ 
tation impact on surface formation. 



Zusammenfassung 


Wie beeinflusst die Vegetationsverteilung die Erosion und Oberflachenfor- 
mung in den verschiedenen Okozonen Chiles? Um diese Frage zu beantwor- 
ten, ist es notwendig grundlegende Kenntnis iiber den Lebensraum, das Vor- 
kommen und die Dynamik von Pflanzenarten zu besitzen. Diese Arbeit ver- 
wendete sowohl frei verfiigbare Sentinel und Landsat Satellitenbilder sowie 
das neue TanDEM-X digitale Gelandemodell (DGM) zusammen mit ver¬ 
schiedenen Ansatzen des maschinellen Lemens (Random Forest, Support 
Vector Machines). Die Ziele waren Karten der Landnutzung und Landbede- 
ckung (LULC) innerhalb von vier Untersuchungsgebieten entlang eines kli- 
matischen Gradienten abzuleiten und Vegetation mittels des normalisierten 
differenzierten Vegetationsindex (NDVI) zu beobachten. NDVI Zeitserien 
basierend auf Einzugsgebiet-Mittelwerten in einem Mediterranen (La Cam- 
pana Nationalpark) und einem humid-temperierten (Nahuelbuta National- 
park) Untersuchungsgebiet wurden statistisch mittels des BFAST („Breaks 
For Additive Seasonal and Trend") Algorithmus untersucht (Verbesselt et al., 
2010) und daraufhin in klassengenaue und hohenspezifische Zeitserien unter- 
teilt. Die Kartierangsergebnisse fielen zugunsten von Random Forest aus, 
welcher die auf Support Vector Machines basierenden Klassifikationen in 
jedem der vier Untersuchungsgebiete iibertraf. Die Gesamtgenauigkeiten der 
fertigen Karten waren stets iiber 85%, wobei eine Flache von 10 000km 2 
abgedeckt wurde und 19 verschiedene LULC Klassen unterschieden wurden. 
Die Zeitserienanalyse stellte eine signifikante Auswirkung des 2015/16 El 
Nino auf das NDVI-Muster in La Campana fest, konnte jedoch keine Belege 
in Nahuelbuta finden. Vegetation in beiden Untersuchungsgebieten zeigte 
eine andauemde Erholung von der zwischen 2011 -2015 stattgefundenen 
zentralchilenischen Mega-Diirre (MD). Die Erkenntnisse konnen zu einem 
besseren Verstandnis der Wirkung des Klimas auf die chilenische Vegetation 
beitragen und eine Grundlage darstellen fur Landschaftsentstehungsmodelle, 
welche den Einfluss der Vegetation auf die Oberflachenformung beriicksich- 
tigen. 



1 Introduction 
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With its elongated shape, Chile is spanning several climatic zones from the 
arid climate of the Atacama Desert in the north to the subpolar oceanic cli¬ 
mate of the Patagonian steppe in the south. Along with this climatic gradient 
goes a spectacular gradient of vegetation including such disparate communi¬ 
ties like the coastal Lomas formation of low growing shrubs and cacti on the 
border of the Atacama and the Valdivian temperate rainforests in southern 
Chile, which are home to the impressive but highly endangered monkey 
puzzle trees {Araucaria Araucana). The EarthShape project as part of which 
this thesis was prepared, focusses on four core study sites, which are located 
inside four different eco zones of Chile. It is the project's pivotal research 
question to study the relationship between vegetation and surface formation 
across this climatic gradient and challenge the common geoscience paradigm 
that holds that landscape evolution is conditioned solely by climate and to¬ 
pography. To do so it is inevitable to have a fundamental knowledge of plant 
and ecosystem distribution but also on the dynamics of vegetation over time. 
Changing climate and its consequences for the Chilean environments add 
further relevance to this necessity. 

At this stage it is widely accepted that Chilean ecosystems are on a constant 
threat due to changing climate conditions (IEB et al., 2010). The urgency of 
already existing and upcoming issues associated with the changes in the 
Chilean environment are documented by adaptation and action plans au¬ 
thored by the Ministry of the Environment (MMA, 2014; MMA 2017). Pre¬ 
vious studies have made major advances in understanding processes of native 
forest loss, landscape fragmentation and the relationship between primary 
production and rainfall in the Chilean environment (La Maza et al., 2009; 
Miranda et al., 2015; Heilmayr et al., 2016; Schultz et al., 2016). However, 
as Parmesan (2006) noted in her study, only comparably less biological stud¬ 
ies come from Chile and South America in general, leaving a significant 
knowledge gap for a continent that has fundamental importance for global 
biodiversity, comprising some of the countries with highest species diversity 
worldwide. The current state of large-scale to nation-wide vegetation and 
land cover maps in Chile shows a lack of studies that took advantage of the 
latest generation of satellite missions. Existing land cover and land use 
(LULC) maps have variable scale, while none exceeds 30 m mapping resolu¬ 
tion (Schulz et al., 2010; Altamirano et al., 2013; Zhao et al., 2016) and in- 


© Springer Fachmedien Wiesbaden GmbH, part of Springer Nature 2019 
C. J. Bodinger, Remote Sensing of Vegetation, BestMasters, 
https://doi.org/10.1007/978-3-658-25120-8_l 




2 


1 Introduction 


vestigates the changes in Chilean vegetation using dense satellite imagery 
time series. Both is now possible with the new Copernicus mission's Sentinel 
satellites and the publicly available long-term archive of the Landsat pro¬ 
gram, which together allow for large-scale mapping and monitoring in un¬ 
precedented resolution. 

This study compliments previous work by providing LULC maps covering 
~10 000 km 2 in transects from the coast up to the Andean mountains includ¬ 
ing each of the four EarthShape primary research areas. These maps have a 
special focus on vegetation to supplement the EarthShape projects focal ob¬ 
jectives and help to answer the question if and how vegetation is changing 
with altitude, latitude and climate in the Chilean Coastal Cordillera. The 
construction of these maps was done using state-of-the-art machine learning 
approaches and a comprehensive set of predictors (multi-temporal optical, 
radar and TanDEM-X DEM derived features). The findings helped to test the 
hypothesis that (high-resolution) DEM and multi-temporal features are (1) 
valuable auxiliaries for vegetation mapping and (2) whether S-l radar fea¬ 
tures can significantly contribute in this classification framework. Time se¬ 
ries of monthly NDVI was produced and analyzed statistically using the 
BFAST (Breaks For Additive Seasonal and Trend) algorithm (Verbesselt et 
al., 2010). This allowed for a detailed consideration of seasonal and multi¬ 
year trends as well as disturbing events that alter the pattern of vegetation 
over time. This study further hypothesized that climate phenomena such as 
El Nino Southern Oscillation (ENSO) have a strong and immediate effect on 
vegetation growth and constitution in Chile and that changes in these plant 
traits can be monitored by NDVI time series. Based on this more general 
hypothesis, this study further tried to explore if the described relationship 
with climate is varying with latitude but also with vegetation species and 
height. 



2 State of Research 
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The EarthShape project, as part of which this master thesis was written, does 
not have a strong remote sensing focus yet. A short introduction into the 
EarthShape projects framework is given pointing out on the major research 
objectives. Based on these, deficits in present research inside the project 
itself and remote sensing in Chile as a whole are outlined. Finally, objectives 
for this study are formulated facing some of the described research gaps. 


2,1 The EarthShape Program 

EarthShape follows an interdisciplinary research approach connecting differ¬ 
ent scientific disciplines in the Geosciences, Ecology, Soil Sciences, Hydrol¬ 
ogy, Microbiology, and Geography. Scientific studies on landscape evolution 
have a strong focus on climate (modulating erosion) and tectonics (fonning 
topography) (Dixon and Thom, 2005; Gallagher et al., 2008). However, the 
consequences of climate change on a larger time scale than just a few centu¬ 
ries remain largely unresolved. In terms of the effects climatic variation has 
on topography, much research has been spent on the influences of increased 
or intensified rainfall and temperature modulated weathering, always applied 
to different bedrock. Also, topographic controls on biodiversity have been 
widely studied (Coblentz and Riitters, 2004). Less is known about the role 
biosphere has on topographic evolution of a landscape. The principal objec¬ 
tive of EarthShape is to explore the role of biologic processes in landscape 
evolution, such as microorganisms and plants on the formation and stabiliza¬ 
tion of soils. The Research efforts try to explore how these processes shaped 
the surface in the past in order to parameterize these biologic controls on 
topography for the future climate and biologic changes. Emphasis is placed 
on non-anthropogenic geomorphic processes, but human influence cannot be 
fully excluded. The project has a strong collaboration with leading Chilean 
scientists from several disciplines (https://esdynamics.geo.uni-tuebingen.de/ 
earthshape). 
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2 State of Research 


2.2 Previous Work & Research Gaps 

Given the framework of the EarthShape project and the present state of re¬ 
mote sensing science, some deficits in research emerge. The principle objec¬ 
tive of EarthShape is to investigate the role of the biosphere on formation of 
topography. Studies have been carried out with very-high resolution using 
lidar and drone data to investigate plant communities, properties and single 
individuals inside the protected areas (Kleinohl, 2017; Senn, 2017). But for 
research on earth surface dynamics also the surrounding country is of im¬ 
portance as the genesis of the landscape is influenced on a larger scale. How¬ 
ever, data about the state of the environment right outside the borders of the 
protected areas is limited. 

Up-to-date information on the distribution of natural and planted forests is of 
great importance to biodiversity conservation. Chilean ecosystems are on a 
constant threat due to changing climate conditions (IEB et al., 2010). Espe¬ 
cially semi-arid regions are expected to experience significant variation in¬ 
cluding a shorter and delayed greening-up season (La Maza et al., 2009). 
Humans increase the force on natural biomes through clear-cutting of natural 
forests and introduction of exotic species. For some while, Chile was regard¬ 
ed as one of the few countries in Latin America to have reversed persistent 
forest loss. But recent studies do now suggest that the increase in forest cover 
is caused by forest plantations of mainly Pinus radiata and Eucalyptus spe¬ 
cies, which are both not native to Chile. Some regions noticed a conversion 
from natural forest to plantations of as much as 16% (Heilmayr et al., 2016). 
For this reason, plantations together with agricultural expansion are some of 
the major drivers behind the loss of natural forest biomes. Multi-temporal 
analysis was carried out to monitor trends like native forest loss and land¬ 
scape fragmentation in south-central Chile (Miranda et al., 2015; Heilmayr et 
al., 2016; Schultz et al., 2016). Previous studies show the complexity of fu¬ 
ture changes in Chilean ecosystems (Alarcon and Cavieres, 2015), which 
makes time series analyses appear more suitable to monitor small-scale 
changes and trends of ecosystems as a consequence of climate change. How¬ 
ever, in Chile detailed time series analysis as available for parts of the Ama¬ 
zonian rain forest are missing so far (Rufin et al., 2015). 

A land cover map with special focus on vegetation communities is available 
to the public through a WMS interface and downloads provided by the min¬ 
istry for agriculture (MINAGRI) (ide.minagri.gob.cl). The mapping scale is 
regionally variable between 1:30 000 and 1:250 000. Only 4 classes attribut¬ 
ed to vegetation are distinguished: 1. Agricultural fields; 2. Temperate grass- 
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lands, savannas, and shrablands (“Praderas y Matorrales”); 3. Forests; 
4. Wetlands. This classification scheme has its limitations in terms of classes 
such as “natural forest”, which must be defined differently as latitudinal 
distance increases. In the northerly regions of Chile, “natural forest” are 
marked by the occurrence of Prosopis trees (i.e. a drought-resistant species). 
More southern, Mediterranean and temperate forests are composed of a 
greater variety of species (evergreen trees like Araucaria Araucana and de¬ 
ciduous such as Nothofagus). The cadaster is based on information from 
1997-2011 and no major actualizations have been made since then 
(www.conaf.cl/nuestros-bosques/bosques-en-chile/catastro-vegetacional). 

A leading remote sensing study on Chilean land use and land cover was pub¬ 
lished by Zhao et al. (2016). Their maps with full coverage of Chile are 
based on Landsat and Modis data resulting in 30 m spatial resolution (i.e. 
1:60 000 mapping scale). The classification scheme includes different levels 
of detail and classification accuracy of the final maps is rather low, especial¬ 
ly for the higher levels of detail (level 1: 80%; level 2: 73%; level 3: 59%). 
Chile covers a whole climatic gradient resulting in a heterogeneous land¬ 
scape, with low vegetation cover in semi-arid to arid conditions. So, it is 
obvious that mixed pixels in 30 m resolution images might become an issue 
in those regions, with the vegetation signal being suppressed by bare soil and 
bedrock. Other maps such as the European Commission’s vegetation and 
land cover maps for South America are even coarser and older (Sul et al., 
2002; Eva et al., 2004) or are spatially restricted with a limited number of 
classes (Schulz et al., 2010; Altamirano et al., 2013). 


2.3 Study Objectives 


Out of the research presented in chapter 2.2, 3 focal hypotheses have 

been derived: 

1. The vegetation is distributed heterogeneously over the different eco 
zones in Chile along a climatic North-South gradient, but also along a 
West-East gradient dependent on the altitudinal belts and distance 
from the coast. 

2. The ENSO phenomenon has a significant impact on the vegetation, but 
the influence is uneven depending on latitude, species and altitude. 

3. DEM features are valuable auxiliaries in vegetation mapping. The 
high-resolution TanDEM-X DEM can compete with established global 
DEMs. 
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2 State of Research 


To investigate hypothesis 1, a vegetation-focused LULC classification is 
produced for a total area of 9503.4 km 2 . Two state-of-the-art machine learn¬ 
ing algorithms (Random Forest, Support Vector Machines) are examined in 
order to investigate, which classifier is superior when high-dimensional da¬ 
tasets from multiple sources (optical, radar, topographic) meet with limited 
training data. Study sites are distributed along a latitudinal and climatic gra¬ 
dient that covers a majority of the main vegetation communities in Chile. 
Every study area comprises also an altitudinal gradient ranging from the 
coast to as far as 150 km into the continent up to the Andes (max. elevation: 
5970 m). Fleight specific analysis is expected to give new insights into state 
and dynamics of vegetation in Chile with implications for climate change. 
The resolution of the maps (10x10 m) with a mapping scale of 1:20 000 is 9 
times higher than those published by Zhao et al. (2016). Apart from classes 
dedicated to vegetation also other (anthropogenic) land use is assessed and 
quantified. However, it is beyond the scope of this study to discuss these in 
great detail. Also, no comprehensive discussion of Chilean plant species is 
given. The reader is referred to Moreira-Munoz (2011) for a detailed over¬ 
view. Species and communities described herein are relevant to the mapping 
scale (i.e. pixel size of the Sentinels). 

To approach hypothesis 2, the impact of climatic and natural forces on plants 
is investigated with a time series analysis. In the frame of this master thesis it 
is impossible to compile land cover maps for multiple points in time. There¬ 
fore, monthly NDVI (Normalized Difference Vegetation Index) as an index 
of plant constitution is calculated inside the most dynamic study sites (La 
Campana, Nahuelbuta), which are exposed to climate phenomena (ENSO) 
and other disturbances (grazing of cows). To investigate if these impacts and 
the resulting pattern of NDVI over time is varying with factors like height 
and species, the time series are split into class-wise and altitude specific time 
series based on the DEM and the classification map produced before. The 
results are analyzed statistically to distinguish between seasonal and long¬ 
term trends as well as breaks in the data that point out on disturbing events 
(Verbesselt et al., 2010). 

Special effort is put into the pre-processing and classification process to ad¬ 
dress hypothesis 3 and achieve a high classification accuracy. Following the 
tasks formulated by Dingle Robertson et al. (2015) for future radar+optical 
LULC remote sensing, variable importance (VI) measures from Random 
Forest (Breiman, 2001) provide insights into the usefulness of certain varia¬ 
bles for LULC mapping. VI is investigated for the relative contribution of 
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predictor variables such as elevation in each study site and per class. The 
newly available TanDEM-X DEM with 12 m spatial resolution is used and 
expected to increase quality of results of both pre-processing and classifica¬ 
tion compared to lower resolution SRTM (Shuttle Radar Topography Mis¬ 
sion) DEM. Different topographic correction methods are applied to S-2 data 
and compared with those of the Sen2Cor algorithm, which was developed 
solely for S-2 pre-processing and whose products are used in the classifica¬ 
tion process. To assess the benefit obtained from topographic correction the 
best performing classification is compared to the best performing classifica¬ 
tion with only atmospheric correction applied. An independent set of valida¬ 
tion data in the difference zone between both classifications is finally used to 
quantify the difference in accuracy following a approach presented by 
Vanonckelen et al. (2013). 
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Figure 1 : Overview on the part of Chile relevant to this study (1.). The four extend¬ 
ed study areas are shown in darkgreen above the mosaiced 12 TanDEM- 
X DEM tiles. Extent of Chile is shaded and the location of the capital 
Santiago de Chile is given by a star. On the right, altitude transects 
through all four study areas are shown. Altitude on the y-axis is in 
m.a.s.l. and the x-axis represents km from the coast. Rectangles mark the 
position of the catchment 

The research efforts focus around four primary areas where all studies work, 
and two secondary areas where a subset of projects are working. This study 
focusses on the primary areas which are presented in the following subsec¬ 
tions. The study sites have been chosen based on a set of criteria making 
them especially valuable for quantifying biotic interactions with topography: 
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3 Study Areas 


1. All of them are situated in the Chilean Coastal Cordillera between 
25° and 40°S. 

2. They have a similar kind of bedrock and glaciation has never 
reached them even during the Last Glacial Maximum (LGM). 

3. Anthropogenic influence can theoretically be neglected as all study 
areas are either long-since natural reserves or national parks. 

Note that the catchments used in this study differ from the outlines of the 
protected areas and are defined on a hydrological basis. From north to south 
a strong climate gradient characterizes coastal Chile. Precipitation gradually 
increases while temperature decreases the further south a study area is locat¬ 
ed. Following the climate gradient, a vegetation gradient goes along, with 
nearly zero percent vegetation cover in Pan de Azucar up to one hundred 
percent canopy cover and dense humid forests in Nahuelbuta. Flowever, this 
study goes far beyond the borders of the core study sites. The reasons are 
mentioned in previous chapters (s. chapters 2.2 and 2.3) and include the 
monitoring of human impact as well as a more comprehensive approach on 
vegetation mapping. As far as the DEMs allow, study areas are stripes cross¬ 
ing from the coast to the border of the Andes with up to 150 km length and a 
common width of 10 km. In the following descriptions a focus is laid on the 
focal EarthShape catchments. 


North Snow limit Upper vegetation limit South 
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Figure 2: Distribution of vegetation formations along the altitudinal profile accord¬ 
ing to Schmithiisen (1956). Rectangles mark the approx, position of the 
study areas (from left to right): 1 st = Pan De Azucar; 2 nd = Santa Gracia; 
3 rd = La Campana; 4 th = Nahuelbuta. Modified after Moreira-Munoz 
( 2011 ). 
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3.1 Pan De Azucar 

Parque Nacional Pan de Azucar (center 26°08’26” S 70°30’52”; catchment 
area: 22.99 km 2 ) is the northernmost study area, located ~10km away from 
the pacific coast, on average ~600 m a.s.l. (Figure 1). For Antofagasta annual 
rainfall is 2 mm and mean temperature 16.5°C. Large seasonal discrepancies 
are missing. The coldest month is July with 13.4°C and the hottest January 
(20.2°C, austral summer). The park is situated within the Tropical macrobio¬ 
climate and the desert core/desert scrub formations according to Schmithusen 
(1956) as part of the Atacama desert (Figure 2). Different species of cacti, 
shrub and stands of planted Prosopis trees (seldom natural) exist close to the 
coast. Plants not watered by men are highly dependent on humidity support 



70°32'24.00"W 70°30'0.00"W 


Figure 3: The Pan De Azucar catchment (upper map) and the study area (small 
lower map). 
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through coastal fog (Rundel and Mahu, 1976). In the park itself, trees are 
missing, sparse shrubs and c acti preferentially occur on the valley floor. In 
the center of the Atacama further east, vegetation is almost completely lack¬ 
ing (Moreira-Munoz, 2011). 

Chile uses the USDA soil taxonomy. To avoid confusion with WRB and 
other taxonomic reference systems, in the following the soil description en¬ 
sues based on their attributes rather than their taxonomic name of soil order. 
Soils in Pan De Azucar are in an early stage of development or often com¬ 
pletely missing. Where sediment accumulates, skeletal and initial soils can be 
found. Salty soils form when ground water is ascending and accumulating 
salts close to the surface. These soils also form in adjacency to salt flats 
(Casanova et al., 2013) such as the smaller salty lake -55 km away from the 
coast, which is partly included in the study area. In the catchment, Triassic 
granitoides form the bedrock with sediments of Holocene and Pleistocene 
age in areas of accumulation. Metasedimentary Paleozoic rocks locate closer 
to the coast and younger Mesozoic volcanic rocks further inside the desert 
(Moreno and Gibbons, 2007). 


3.2 Santa Gracia 

The Reserva Natural Santa Gracia en La Serena (center 29°44’09” S 
71°08’03” N; catchment area: 82.3 km 2 ) is situated on the transition between 
Tropical macrobioclimate and Mediterranean-type macrobioclimate. In con¬ 
trast to Pan De Azucar, the study area reaches far into the Andes above 5000 
m a.s.l. (Figure 1). Mean annual rainfall is 95.9 mm (1919-2013) and mean 
temperature is lower than in Pan De Azucar (13.6°C). Winter rainfall (June- 
August) accounts for 70% of the total annual rainfall. Since the beginning of 
measurements, annual rainfall has decreased by about 6 mm. This puts the 
region on the risk of desertification, while flooding after intense rainfall 
events caused by El Nino is another threat to flora and fauna (Soto et al., 
2017).Vegetation formations comprise xerophytic scrubs, replaced by desert 
scrubs above -1000 m and low Andean vegetation above -2000 m (Figure 
2). 

The area is famous for the coastal so-called “Lomas” formation, a rich as¬ 
semblage of shrub species ( Eulychnia , Nolana, Heliotropium, Tetragonia, 
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Figure 4: The Santa Gracia catchment (upper map) and the study area (small lower 
map). 


and Euphorbia ) highly dependent on fog and to some extent El Nino (Rundel 
and Dillon, 1998). During El Nino events the desert experiences a sudden 
grow of flowering plants. This rare phenomenon is known as the “blooming 
desert” (span, “desierto florido”), which occurred in August 2017 after in¬ 
tense rainfall earlier in the year. With a close look on Figure 4, shmbs are 
visible in accumulation areas and valleys incised by highly erosive but rare 
rainfall. Planted and naturally grown Prosopis trees gather together with 
cacti and scrub species. These communities are visible even in 10 m resolu¬ 
tion imagery with a clearly discernable vegetation signal. Irrigated agricul¬ 
ture can be found surrounding La Serena City but also far away from the 
coast along and close by the route 41 crossing Chapilca and Guanta. Intense 
mining activities have shaped this part of Chile for centuries. Immediately 
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west of the reserve lies the iron mine of El Romeral, a depression easily visi¬ 
ble in the DEM. Bedrock is formed by large plutonic rocks and younger 
sediments. Further east, volcanic rocks of younger Mesozoic origin like bas¬ 
alt and andesite dominate (Moreno and Gibbons, 2007). 

Soil formation is slightly more advanced compared to the center of the Ata¬ 
cama further north. Coarse textured soils exist on dunes close to the coast or 
in colluvial sediments. Soils of granitic origin have a pH between 6.3 and 8.0 
and surface organic matter content <2% and are assumed to be resistant 
against erosion (Casanova et al., 2013). Andosols, especially on steep slopes, 
are more prone to erosion. 


3.3 La Campana 

La Campana (center 32°56’59” S 71°03’53” N; area: 49.1 km 2 , i.e. the Ocoa 
catchment) was declared a national park in 1966 and Biosphere Reserve in 
February of 1985 jointly with Lago Penuelas National Reserve (Rundel and 
Weisser, 1975; Hauenstein, 2012). On his visit to Chile in 1834, Charles 
Darwin climbed the Cerro La Campana an 1880 m high peak in the parks 
west (Figure 5). The climate belongs to the temperate-mediterranean with 
strong seasonality. The average annual temperature is 13.5°C; the average 
maximum is 17.1°C and the average minimum is 9.4°C. Especially in higher 
elevations, temperature might occasionally fall below 0°C. Annual precipita¬ 
tions is 656 mm (Hauenstein, 2012). According to (Manriquez, 2002), cli¬ 
mate of the park is classified into three seasons: 1. October to March: period 
of dry and hot weather, 2. April to August: rainy period with unstable climate 
conditions and 3. September: spring period with a lot of fog and cyclonic 
disturbances. Theoretically, formations of sclerophyllous matorral, xerophyt- 
ic and desert scrubs dominate the landscape (Figure 2). However, La Cam¬ 
pana is home to probably the highest diversity of vegetation communities 
across all study areas. Hauenstein (2012) lists 420 species of vascular plants, 
of which 49% are of native origin, 28% corresp ond to endemic species and 
23% to introduced species. All major biotic communities of central Chile are 
present in the park: Nothofagus forest, hygrophilous forest, sclerophyll for¬ 
est, matorral, bamboo thicket, succulent scrub, high altitude communities, 
and palm forest (Rundel and Weisser, 1975). Many species are endemic to 
Chile such as the Chilean palm (Jubaea chilensis), which grows in the valley 
of Ocoa. Centered around 400 m, hygrophilous forests develop in the valleys 
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Figure 5: The Ocoa catchment (upper map) and the study area (small lower map). 


above Olmue. These forests are strongly influenced by coastal fog and show 
a rich diversity of different woody species (e.g. Drimys winteri (canelo, holy 
tree of the Araucanian Indians), Crinodendron patagua (patagua), Persea 
lingue (lingue) and the climber Proustia pyrifolia (parrilla blanca)) and 
mosses and lichens (Rundel and Weisser, 1975).Vegetation cover varies with 
biotic communities. Nearly 100% are reached in the evergreen sclerophyll 
forests of Peumo and Boldo types. This community shows a strong prefer¬ 
ence of southern, south-eastern and eastern aspects. In general, higher plants 
preferably grow on south-facing aspects rather than on north-facing. South¬ 
facing slopes between 1100 and 1500 m on El Roble, a 2220 m high peak in 
the parks east, also host the northermnost stands of Nothofagus (a deciduous 
oak) trees. North-facing slopes show a rich array of annual species and char¬ 
acteristic bromeliads, together with the cactus Trichocereus chiloensis 
(Moreira-Munoz, 2011). Coverage at these sites and on the most exposed and 
flat areas where Prosopis and acacia caven (bromeliads) grow, is considera¬ 
bly lower. Matorral dominated by Lithraea caustica (litre), Quillaja sapona- 
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ria (quillay) has variable cover between 50-100%. This community is an¬ 
other typical biotic formation in the La Campana region, with low, evergreen 
shmbs between 1 and 2 m height (Rundel and Weisser, 1975). Environmen¬ 
tal degradation threatens these ecosystems. In the park grazing with cows, 
although officially forbidden, causes considerable damage. Plant composi¬ 
tions outside the park is heavily influenced by anthropogenic land use. Be¬ 
sides agriculture, orchards (i.a. olives, oranges, mandarins) and vineyards, 
plantations of exotic species are common (coniferous Finns radiata and ev¬ 
ergreen Eucalyptus species). Most often plantations are characterized as 
regular stands of individuals of different age. 

The bedrock in the catchment is once more formed by igneous intrusive 
rocks (Diorites, Monzodiorites, Granodiorites) and younger sediments. Out¬ 
side the catchment, geology is of similar complexity than in the other study 
sites. Volcanic rocks of different ages occur especially towards the Andes, 
while Pleostocene glacifluvial and Holocene alluvial deposits are common in 
accumulation areas (Moreno and Gibbons, 2007). According to Bonilla and 
Johnson (2012) soils in this region are rather shallow and slightly acidic. 
Nonetheless, fertility of soils formed out of colluvial material is usually suf¬ 
ficient for agriculture. Towards the Andes with increasing elevation, soils are 
similar to the shallow soils of Santa Gracia and desertification becomes more 
important (Casanova et al., 2013). 


3.4 Nahuelbuta 

Parque Nacional Nahuelbuta (center 37°48’16” S 72°59’02” N; catchment 
area: 28.7 km 2 ) belongs to the temperate macrobioclimate with Mediterrane¬ 
an influence (Figure 6). The park was founded in 1939 and is a reserve for 
the endangered coniferous Monkey puzzle tree (Araucania Araucana ). The 
term “Nahuelbuta” stems from the Mapuche indian people and can be trans¬ 
lated with “big tiger”, presumably named after the mountain lion living in the 
park. Average temperatures range between -1°C in winter and +9°C in sum¬ 
mer. Annual precipitation is between 1500-2500 mm, occasionally up to 
5000 mm (Zamorano-Elgueta et al., 2012). 

The characteristic pure Araucania forests are limited to the highest points in 
the national park where bedrock is outcropping (Endlicher and Mackel, 
1985). Mixtures of Nothofagus dombevi (Coihue) and Araucania occur on 
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Figure 6: The Nahuelbuta catchment (upper map) and the study area (small lower 
map). 


less rocky summits. Robles grow on the warm slopes, while the humid 
gorges (span, “quebrada”) are dominated by Nothofagus dombeyi. Due to 
legal and illegal harvesting, almost no Araucania trees exist outside the park. 
Native forests in the Nahuelbuta range have suffered great losses following 
the expansion of commercial plantations, agricultural activities and urban 
and industrial sprawl. As studies suggest, although being protected the trees 
are still threatened by constant direct and indirect human disturbances (fires, 
cattles, selective logging of forest species associated with the conifer) (Za- 
morano-Elgueta et al., 2012). 

Using historical documents, archaeological evidence and environmental 
modelling, Lara et al. (2012) show a decrease of native forests from 11.3 
million ha ca. 1550 at the Spanish conqueror arrival, to 5.7 million ha in 
2007 (50.4% of the original area). Outside the Nahuelbuta range natural for¬ 
ests have survived especially at steep slopes. This is the case north of lake 
Lanalhue on south-exposed slopes, which are in general steeper than those 
facing to the north. Pasture and agricultural fields dominate the easier acces- 
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sible north-facing slopes. West of Canete (Rucaraqui area) the landscape is in 
transition to the littoral plains, characterized by typical box-shaped valleys. 
On the valley heads and around seasonal pans, natural vegetation is com¬ 
posed of a forest rich in Chilean hazelnut ( Gevuina avellane, called Avellano 
by inhabitants) and Myrtaceae. Steeper slopes are home to Peumus boldus 
and Aextoxicon punctatum (Olivillo) forests with dense understory. River 
banks are covered with dense bamboo ( Chusquea quila) thicket, liana and 
low evergreen or deciduous shmbs (Ugni molinae, Senecio yegua, Jovelana 
violacea). The valley floor is often used for grazing or agriculture (wheat, 
potatoes). The up to 6 km wide littoral plains are characterized by dunes and 
woodlands formed by Peumus boldus. Removal of the natural vegetation 
leads to reactivation of formerly stabilized dunes. To address this issue, CO- 
NAF (National Forest Corporation) has decided in 1985 that dunes should be 
planted with grasses of the genus Ammophila arenaria and Pinus radiata 
plantations. Both species are tolerable to sandy soils (Endlicher and Mackel, 
1985; Moreira-Munoz, 2011). 

The Nahuelbuta range is formed by intrusive (mainly granitic, granodiorotic 
and dioritic) and metamorphic rocks. Fluvial and coastal sandy sediments 
occur closer to the sea and on the littoral plains. Rocks of volcanic origin 
become frequent towards the Andes (Moreno and Gibbons, 2007). Soils in 
the park originating from metamorphic material are thin to moderately deep 
(15-180 cm). The texture is a sandy clay loam with strong to moderate acidi¬ 
ty (pH 3.0-5.5) and low nutrient level. Erodibility is high as can be seen out¬ 
side the park where gullies have formed after removal of the forest cover 
(Zamorano-Elgueta et al., 2012). 
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The following subsections shall provide an overview on the data used in this 
study as well as its sources. Furthermore, the techniques applied to this pri¬ 
mary data to derive secondary data and the methods for classification and 
assessing the accuracy are presented in some detail. 


4.1 Data 

This study is based upon satellite data and derived products (DEM) as well 
as (vector) datasets accessed through the internet (s. subsection 4.1.3). All 
data, except from TanDEM-X high-resolution DEM, are open source and can 
be downloaded without restrictions. Please note that some datasets (e.g. Sen¬ 
tinel and Landsat data as well as data provided by Global Biodiversity In¬ 
formation Facility (GBIF)) require preliminary registration. 

4.1.1 Sentinel Optical and Radar Data 

The LULC classifications are based upon imagery from the Sentinel satellites 
of the European Union's Copernicus program (fonnerly Global Monitoring 
for Environment and Security (GMES)). S-l is a two-satellite constellation 
since the launch of S-1B on the 25 th April 2016. These C-band SAR (Syn¬ 
thetic Aperture Radar) satellites provide imagery in either HH+FIV, VFI+VV, 
HH or W polarization covering the entire earth every six days with the spa¬ 
tial resolution depending on the acquisition mode and the level of processing. 

In this study Interferometric Wide-swath (IW) Level-1 Ground Range De¬ 
tected (GRD) imagery of ~10m spatial resolution are used. GRD data is fo¬ 
cused SAR data that has been detected, multi-looked and projected to ground 
range using an Earth ellipsoid model (ESA, 2013). Polarizations are dual- 
polarized VH+VV. SAR data provides information about surface roughness, 
soil moisture, object shape, orientation and texture. 

S-2 is an optical two-satellite constellation with the launch of the second 
satellite on 7 th March 2017. The S-2 satellites have 290 km swath width re¬ 
sulting in a 10-day revisit cycle and near global coverage (56° South to 83° 
North). The 13 bands cover a spectral range from 443 — 2190 nm with three 
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bands in the vegetation-sensitive red edge region centered at 705, 740 and 
783 nm. The 10, 20 and 60m resolution bands and a sun-synchronous orbit at 
786 km altitude with a descending node allow for continuity with SPOT and 
Landsat. A local overpass time of 10:30 a.m. also minimizes cloud cover and 
ensures suitable sun illumination (ESA, 2010). S-2 imagery is provided in 
Top Of Atmosphere (TOA) Level-lC geo-coded reflectance with a sub-pixel 
multi-spectral and multi-date registration. Level-2A Bottom Of Atmosphere 
(BOA) reflectance is not yet distributed but can be derived using ESA's 
Sen2Cor processor or other image processing software. Sen2Cor is available 
as a third-party plugin to the S-2 toolbox or ESA's SNAP software (Louis et 
al., 2016; Gascon et al., 2017). 

The Sentinel missions provide a unique archive outperforming the Landsat 
program at least in terms of resolution. Studies predicted that the Sentinels 
will lead to a better understanding of bio-geophysical variables and thus 
fluxes and interactions between the surface and the atmosphere once longer 
dense time series are available (Berger et al., 2012; Sakowska et al., 2016). 
All Sentinel data used herein were downloaded from the official data hub 
hosted by ESA (www-scihub.copernicus.eu/dhus). Imagery for classification 
maps was acquired in the summer months, which is the most representative 
time in the year regarding vegetation. Images for the seasonal NDVI, which 
were used as a predictor in classification, were acquired once in the begin¬ 
ning of every season. Cloud cover is low across all optical images, only clas¬ 
sified images in Pan de Azucar and Santa Gracia have little cloud cover 
(>10%). Coverage of the optical scenes was except for Pan De Azucar insuf¬ 
ficient for the whole study areas. Therefore, two images were acquired and 
mosaiced. Due to large discrepancies in grey values and illumination condi¬ 
tions between the two scenes in La Campana, those were processed inde¬ 
pendently resulting in two classifications. In the following La Campana 1 
denotes derived products and the scene closer to the coast, while La Campa¬ 
na 2 refers to the scene covering the inland part. Acquisition of radar images 
ensued with smallest possible time gap to S-2 scenes. Only in Santa Gracia 
two radar images were necessary to cover the whole study area. 

Both radar and optical data are used in combination with topographic attrib¬ 
utes to improve mapping accuracy. Combinatory usage and data fusion have 
improved results in various studies and different environments (Dingle Rob¬ 
ertson et al., 2015; Gianinetto et al., 2015). Optical imagery is more often 
used for mapping purposes and has several well-known advantages (e.g. 



4.1 Data 


21 


Table 1: Acquired datasets from S-l and 2. S-2 acquisitions in bold have been used in 
classification, all others are used for seasonal NDVI. In brackets the number 
of tiles necessary to cover the study area. 


Study Area Sensor Acquisition 

Study Area Sensor Acquisition 

Pmi Dc §2 19.01.2017 

Azucar 

20.03.2017 

13.07.2017 

16.09.2017 

SI 15.01.2017 

La 

S2 19.01.2017 (2x) 

Campana 

10.03.2017 (2x) 

18.06.2017 (2x) 

16.09.2017 (2x) 

SI 10.01.2017 

) >anta S2 19.01.2017 (2x) 

Gracia 

09.04.2017 (2x) 

08.06.2017 (2x) 

06.10.2017 (2x) 

SI 15.01.2017 (2x) 

15.01.2017 

Nahuelbuta S2 19.01.2017 (2x) 

20.03.2017 (2x) 

18.07.2017 (2x) 

16.09.2017 (2x) 

SI 15.01.2017 


easier interpretability; certain wavelengths corresponding to vegetation con¬ 
ditions) (Joshi et al., 2016). Unlike optical products, radar sensors are un¬ 
affected by atmospheric conditions apart from short wavelengths and dense 
clouds, penetrate vegetation cover providing information about the under¬ 
story and are able to operate day and night thus increasing temporal availa¬ 
bility. 

4.1.2 TanDEM-X DEM 

TanDEM-X is a radar interferometry mission lead by the German Aerospace 
Center (DLR) with the primary objective to generate a consistent global digi¬ 
tal elevation model (DEM) (Krieger et al., 2007). In this study 12 high- 
resolution TanDEM-X tiles were used for which the application has been 
done during an open call in late 2016. These tiles cover a total land surface 
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area of 91 379 km 2 , which were clipped to the size of the study areas (area: 
10 685.25 km 2 , catchments: 312.82 km 2 ). The overall good performance and 
vertical accuracy of the DEM is certified by various (preliminary) compara¬ 
tive studies indicating an error below 2 m or even less for moderate terrain 
(Rizzoli et al., 2012; Wessel et al., 2014; Balzter et al., 2016). All rasters 
with a nominal spatial resolution of 12 m were resampled to 10 m using bi¬ 
linear interpolation and projected to UTM coordinates to match the spatial 
resolution of the Sentinels. Gaps in the original DEMs (Santa Gracia, espe¬ 
cially La Campana) due to radar shadowing mountainous terrain were filled 
using Inverse Distance Weighting in ESRI ArcMap 10.5. Other DEMs used 
in this study and compared to the TanDEM-X DEM include 30 m resolution 
SRTM DEM tiles and a 2 m resolution Lidar DEM. The Lidar DEM is avail¬ 
able in La Campana only and covers solely the Ocoa catchment. It has been 
generated within the first EarthShape project phase. 

4.1.3 Training & Validation Data 

Training data were compiled from various sources to obtain a database with 
vegetation observations distributed across the whole study area. Sample loca¬ 
tions were extracted from SINIA (Sistema Nacional de Information Ambien- 
tal; www.sinia.cl), the BIEN database (Botanical Information and Ecology 
Network; www.bien.nceas.ucsb.edu/bien) and herbaria listed in the GBIF 
database (Global Biodiversity Information Facility; www.gbif.org). In total 
232 plant locations were collected this way and data older than 2010 was 
rejected. As most of the data holds information about trees, it was decided 
that a time span of 7 years could be accepted. Flowever, a large number of 
data points is located within the La Campana study area (210) and only 12 
are situated in Nahuelbuta, 8 in Santa Gracia and 2 in Pan De Azucar. All 
sample locations were visually checked for validity using high-resolution 
imagery in GoogleEarth. 

Information on non-vegetation classes such as Built-up areas and inland 
water bodies as well as forest plantations and other agricultural sites (e.g. 
orchards, vineyards) is extracted from the geospatial information service of 
the MINAGRI (Ministerio de Agriculture; http://ide.minagri.gob.cl/ 
geoweb/). The service also provides a WMS (Web Mapping Service) that can 
assist during on-screen digitization (http://ide.minagri.gob.cl/maps/Publico/ 
wms?). 
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Based on the point and polygon information of the aforementioned sources, 
polygons were digitized in homogenous areas using expert knowledge and 
high-resolution imagery available through Google Earth. The size of the 
polygons highly depends on the respective class. While water bodies are 
available with huge surfaces in all 4 study sites, especially natural vegetation 
classes are usually fragmented within small patches of only a few pixels size. 
Therefore, additional points were set after the polygons were converted to 
raster and subsequently to point format. This explains the huge discrepancies 
in the number of samples available per class. 

Validation data was collected independently from training data. As with 
training data, sample points are digitized close to known occurrences. All 
data was collected in a way that ensures validation and training data do not 
overlap in order to minimize autocorrelation (Congalton and Green, 2009). 
To cover the spectral range of single classes, data was gathered systematical¬ 
ly over the total image following the approach presented by Vanonckelen et 
al. (2013). Spreading the sampling across the whole area was sometimes 
impossible, because certain classes occur only in specific regions with lim¬ 
ited extent. E.g. the Araucaria Auracana trees in Nahuelbuta, which can only 
be found inside or in immediate proximity to the national park. Due to the 
polygon- and non-randomized sampling it is obvious that spatial autocorrela¬ 
tion could not be fully excluded. Considering the huge area covered in this 
solely remote-sensing based study and the relatively sparse data base, it is 
impossible to avoid this drawback. A detailed overview on the number of 
point samples per class split in training and validation is given in Table 2. 
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Table 2: Number of samples per class and study area split in Training/Validation. 


Class 

# samples 

Class 

# samples 

Class 

# samples 

Pan De Azucar 

La Campana #1 

8 

103/70 

9 

266/139 

1 

296/72 

9 

2219/179 

11 

120/52 

2 

59/100 

11 

144/45 

17 

25/18 

3 

166/93 

20 

98/47 

20 

50/56 

4 

160/78 

21 

1481/44 

21 

67/58 

5 

94/44 

22 

124/58 

23 

48/50 

6 

187/53 


Nahuelbuta 

Santa Gracia 

7 

249/54 

2 

359/271 

9 

181/ 155 

8 

324/71 

5 

227/74 

11 

90/55 

9 

128/45 

6 

677/206 

16 

116/88 

11 

291/63 

7 

512/138 

17 

78/53 

20 

29/39 

8 

945/105 

18 

118/66 

21 

1409/57 

10 

181/47 

19 

60/80 

La Campana #2 

11 

409/105 

20 

55/42 

1 

205/80 

12 

140/100 

21 

140/62 

2 

117/95 

13 

339/59 

22 

81/43 

3 

195/83 

14 

170/73 

23 

58/51 

4 

226/72 

15 

654/132 



5 

72/38 

21 

967/60 



7 

141/122 




4.1.4 Landsat Surface Reflectance 

Landsat Surface Reflectance Level-2 Science Data Products were obtained 
from the USGS (U.S. Geological Survey). Surface Reflectance products were 
corrected for atmospheric influences using a different processing scheme 
depending on the respective Landsat mission. Processing for Landsat 4-7 
using the Landsat Ecosystem Disturbance Adaptive Processing System (LE- 
DAPS) is described in Masek et al. (2006). A modification of this algorithm 
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(Landsat Surface Reflectance Code (LaSRC)) taking advantage of the Land- 
sat 8's Aerosol band is described in Vermote et al. (2016). Atmospheric cor¬ 
rection is mandatory for multi-temporal change detection analysis, especially 
for the kind of time series analysis applied in this study (cp. chapter 4.2.8) 
(Lu et al., 2014). Correction for topographic effects was omitted as vegeta¬ 
tion in the studied Ocoa catchment is distributed mainly in flat to moderate 
slopes. Steep slopes inside the catchment on e.g. Cerro La Campana are 
largely devoid of higher vegetation. More than 210 Landsat NDVI scenes for 
La Campana and Nahuelbuta were ordered from the USGS. Only scenes 
without cloud cover in the catchment were chosen and quicklooks checked 
beforehand. The final number of available NDVI scenes amounts to 207. 
Thereof 160 are located in La Campana (06.08.1997-29.08.2017) and 47 in 
Nahuelbuta (12.04.2017-20.08.2017). The reason for this discrepancy are the 
differences between the NDVI calculated from different Landsat missions 
(see chapter 5.5 for more details). Therefore, it was decided to limit the time 
covered by the Time series analysis (see chapters 4.2.8 and 5.5) to 4 years in 
which Landsat 8 data were available (April 2013 until March 2017). 


4,2 Methods 


Table 3: Overview on the variables used for LULC classification. For Abbrevia¬ 
tions refer to the Abbreviations section of this document. 


Spectral Bands (S-2) 

Blue, Green, Red, Vegetation Red Edge [5-7 + 8A], NIR, 
SWIR [11-12] 

Radar Backscatter 

W, VH 

GLCM 

Mean, Variance, Homogeneity, Contrast, Dissimilarity, 
Entropy, 2nd Moment, Correlation) 

Spectral Indices 

CTVI, GEMI, GNDVI, MSAVI2, NDVI, SLAVI, SR, 

TTVI, IRECI, S2REP, VERRELST 

DEM + derived features 

DEM, Aspect, Slope, LS Factor, Longitudinal Curvature, 
Cross-Sectional Curvature, VDCN, ValleyDepth, TPI, 

TRI, TWI, SWI, Direct Insolation 


The predictor set for land cover classification consists of the spectral bands 
(2-8A, 11+12) and backscatter (VV, VH), 10 GLCM features, 11 spectral 
indices as well as the DEM and 13 predictors derived from the DEM. An 
overview on the predictors is given in Table 3. 
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4.2.1 Sentinel Preprocessing 

The magnitude of the signal received from a satellite sensor is dependent on 
several factors, particularly: reflectance of the target; nature and magnitude 
of the atmospheric interactions; slope and aspect of the ground target area 
relative to the solar azimuth; and angle of view of the sensor and solar eleva¬ 
tion angles (Lantzanakis G. et al., 2017). Atmospheric correction of S-2 im¬ 
agery was performed with the freely available Sen2Cor algorithm, which is 
based on ATCOR (Atmospheric & Topographic Correction) (Louis et al., 
2016; Gascon et al., 2017). Sen2Cor allows for the adjustment of a few pa¬ 
rameters (atmospheric temperature profile and ozone content, Cirrus correc¬ 
tion, Bidirectional Reflectance Distribution Function (BRDF) correction), 
while most of the code is proprietary and cannot be accessed. Additional 
topographic correction is recommended if 5% of the pixels in the region of 
interest have slopes >8% (Main-Knom et al., 2015). Removal of the topo¬ 
graphic effect was performed twice for the January scene in La Campana 1 
using SRTM and TanDEM-X DEM. All other study areas and scenes col¬ 
lected for calculation of NDVI were processed with the TanDEM-X DEM. 
To assess the performance of the new Sen2Cor algorithm in topographic 
correction, 8 traditional methods (Gamma, SCS, Improved Cosine, C- 
correction, C-correction with NDVI, Minnaert, Minnaert with Slope, Cosine) 
were applied and compared based on visual interpretation and descriptive 
statistics. While most of these algorithms are already implemented in R 
(package “Landsat”), a customized script was developed to include the 
NDVI into the C-correction. 

Riano et al. (2003) judge the performance of a correction algorithm by its 
effects on mean and standard deviation of each band. The mean of the cor¬ 
rected image is subtracted from the mean of the uncorrected image. The au¬ 
thors consider a better correction as such, that it maintains the mean in the 
band values. The same applies for the standard deviation (SD). Prerequisite 
to do so is that each image is scaled to the same value range. In this study, 
values were scaled to reflectance between 0 and 1. 

Already focused and projected 10 m GRD data (Ground Range Detected) 
from S-l were pre-processed in SNAP v. 5.0.0 according to the workflow 
shown in Figure 7. After evaluation of preliminary results. Terrain Flattening 
(Small, 2011) was omitted due to the introduction of artifacts in mountainous 
areas. Refined lee filter was chosen for speckle filtering as suggested by (Qi 
et al., 2010; Dingle Robertson et al., 2015). Projection to the specified coor- 



4.2 Methods 


27 



Figure 7: Processing chain for S-l imagery in SNAP. 

dinate system and terrain correction were done with the Range-Doppler Ter¬ 
rain Correction approach in SNAP (Laur H. et al., 2004; Bayanudin and 
Jatmiko, 2016). 

4.2.2 Derivation of Terrain Variables 

Topographic characteristics are important drivers of various regional- and 
local scale geodynamic processes. Therefore, topography-derived parameters 
are frequently used for landscape modelling and classification purposes 
(Balzter et al., 2015; Conrad et al., 2015). 

13 topographic indices were calculated from the TANDEM-X DEM using 
SAGA (System for Automated Geoscientific Analysis) GIS v. 2.3.2 software 
(Conrad et al., 2015)(Table 4). Refer to the SAGA online manual for more 
details to the parameters (www.saga-gis.org/saga_tool_doc/2-2.1). 
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Table 4: Topographic variables calculated from TanDEM-X in SAGA GIS. 


Indices 

Reference 

Indices 

Reference 

Aspect 

(Travis et al., 1975) 

SWI 

(Boehner et al„ 

2002) 

Cross-Sectional 

Curvature 

(Zevenbergen and Thome, 
1987) 

TPI 

(Guisan et al., 1999) 

Direct Insolation 

(Bohner and Antonie, 2009) 

TRI 

(Riley et al., 1999) 

Longitudinal 

(Zevenbergen and Thorne, 

TWI 

(BEVEN and 

Curvature 

1987) 

KIRKBY, 1979) 

LS Factor 

(Moore et al., 1991) 

Valley 

Depth 

Author: O.Conrad (c) 
2012 (unpublished) 

Positive 

Openness 

(Yokoyama et al., 2002) 

VDCN 

(Conrad, 2005) 

Slope 

(Travis et al., 1975) 




4.2.3 GLCM 

Texture measures can provide vital additional information from SAR images. 
The Gray Level Co-occurrence Matrix (GLCM) (Haralick et al., 1973) is 
often used for land mapping purposes as the number of bands a radar system 
is usually more limited than those of optical sensors (Balzter et al., 2015; 
Zakeri et al., 2017). The GLCM matrix computes measures in a defined win¬ 
dow size (here: 5x5 pixels), thereby including the information from neigh¬ 
boring pixels. 10 GLCM measures were computed in ESA SNAP v. 5.0.0, 
they can be ordered in 3 groups: 

• Contrast Group Features: 


Contrast 


Info 1 " 2 fair 

l IWI =n| J 

(1) 


o Dissimilarity 

Tnj =0 Pi,j 1 n ^ij) 


( 2 ) 
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o Homogeneity 

2 / i+(£-y) 2 (^ty) P) 

• Orderliness Group Features: 

o Angular Second Moment (ASM) 

EiEjW,;')} 2 (4) 

o Energy 

VASM (5) 

o Maximum Probability 

Max(Pi,j ) (6) 

o Entropy 

-EiPlog(P(i,;)) (7) 

• Statistics Group Features: 
o Mean 

2 iP *+y (0 ( 8 ) 

o Variance 

SiSy(i — H) 2 P(i,» (9) 

o Correlation 

Hi PxPy (10) 


where p(/, y) is the (/, y)-th entry in a normalized gray-tone spatial depend- 
ence matrix P(i,j)/R; R is the total sum of P; p x (i) = ^ is the i-th 

entry in the marginal probability matrix obtained by summing the rows of 
p(i, j); and p x , p y , o x and o y are the means and standard deviations of p x and 
p y (Zakeri et al., 2017). 
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4.2.4 Spectral Vegetation Indices 

A set of 11 spectral Vegetation Indices (VI) was computed based on S-2 
imagery. These indices comprise 8 traditional indices, which have been test¬ 
ed in various preceding studies. There will be no thorough discussion of 
those Vis as this would push the boundaries of this study. The reader is re¬ 
ferred to the references given in brackets in Table 5. 

In addition to the set of well-known indices, 3 newly formed indices are 
used. These Vis have been developed based on simulated S-2 imagery with 
the purpose to take advantage of the unique properties of this new satellite 
mission (Frampton et al., 2013; Verrelst et al., 2015). 

The two indices introduced by Frampton et al. (2013) utilize the RE bands 5 
and 6 to measure Leaf Canopy Chlorophyll (LCC) and Leaf Area Index 
(LAI) even above the common saturation point. As these are both important 
biophysical indices, IRECI and S2REP are integrated into the feature set to 
assess their effectiveness for vegetation mapping. While the IRECI is a com¬ 
pletely new index, S2REP is a modified version of the MTCI (MERIS Ter¬ 
restrial Chlorophyll Index). Due to the slightly different positioning of the S- 
2 RE band 6 on top of the linear part of the RE slope, S2REP should in theo¬ 
ry provide better results than with MERIS (Medium Resolution Imaging 
Spectrometer). 

The third Sentinel index differs from the aforementioned because it uses the 
two bands situated in the SWIR region instead of the RE bands. It was found 
by assessing more than ten thousand artificially constructed indices regarding 
their capabilities for predicting LAI (Verrelst et al., 2015). Like IRECI and 
S2REP, this new unnamed index (referred to as VERRELST in this study) 
differs from many traditional indices because it includes more than two 
bands. 
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Table 5: Spectral Vegetation Indices with formula and source reference. 


Indices 

Formula 

Reference 

CTVI 

(NDVI + 0.5) 
y/\(NDVI + 0.5)| 

(Perry and 

Lauten- 

schlager, 

1984) 

GEMI 

(((NIR 2 - RED 2 ) * 2 + (NIR * 1.5) + RED * 0.5))^ 

^ (NIR + RED + 0.5) ) * ^ 

((NIR A 2 - RED A 2) * 2 + (NIR * 1.5) + (RED * 0.5)) 

(Pinty and 
Verstraete, 

1992) 


(NIR + RED + 0.5) * 0.25 
(RED - 0.125) 

(1 - RED) 


GNDVI 

(NIR - GREEN) 

(NIR + GREEN) 

(Gitelson and 
Merzlyak, 

1998) 

MS AVI2 

(2 * (NIR + 1) - sqrt((2 * NIR + 1) A 2 - 8 * 

(NIR - RED))) 

(Qi et al., 

1994) 

NDVI 

2 

(NIR - RED) 

(NIR + RED) 

(Rouse et al., 
1974) 

SLA VI 

NIR 

(RED + SWIR2) 

(Lymbumer et 
al, 2000) 

SR 

NIR 

RED 

(Birth and 
McVey, 1968) 

TTVI 

1 [NIR - RED] 
l( NIR + RED) + 0.5 

(Thiam, 1997) 

IRECI 

(NIR - RED) 

RE2 

RE 1 

(Frampton et 
al„ 2013) 

S2REP 

(NIR + RED) — RE1 

705 + 35 * , 

(RE 2 - RE 1) 

(Frampton et 
al, 2013) 

VER- 

RELST 

(GREE - SWIR1 - SWIR2) 

(GREEN + SWIR1 + SWIR2) 

(Verrelst et 
al„ 2015) 
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4.2.5 Land Cover Classes 

In total 19 LULC classes were differentiated in this study (excl. shadow, 
clouds and classes later rejected and merged). The decision process for the 
final classes comprised several steps: 

(1) First, a review of existing studies was conducted and some preliminary 
class descriptions were formulated thereof (Schulz et al., 2010; Altamirano et 
al., 2013; Heilmayr et al., 2016; Zhao et al., 2016). 

(2) In a second step, high-resolution imagery was analyzed for visual separa¬ 
bility from other classes. If separability was too low, class merging or dis¬ 
carding was considered. Otherwise spectral separability in Sentinel bands 
and relationships with other predictors (e.g. topography) was considered. 

(3) Finally, Jeffries-Matusita (J-M) distance measures were calculated in R. 
The J-M distance is a function of separability that directly relates to the 
probability of how good a resultant classification will be (Bruzzone et al., 
1995). If separability between classes is low, the classification accuracy will 
likely be limited, while a combination of those classes will lead to a gain in 
map quality. In this study, J-M distance was calculated as separability based 
on all predictors including topographic and radar variables (cp. chapter 4). 

The resulting land cover classes comprise some species-specific classes such 
as Eucalyptus, Pine and Araucaria. This is either due to their occurrence in 
larger homogenous plantations (Pine and Eucalyptus) or due to their signifi¬ 
cance and dominance in some ecosystems (Araucaria in Nahuelbuta National 
Park). In some cases, a detailed separation was impossible, e.g. in the case of 
Mediterranean forests. It was not feasible to sufficiently distinguish sub¬ 
communities like the palm “forests” of La Campana or different kinds of 
shmb communities. 

After classification, some of the classes mentioned above were merged to¬ 
gether. This is due to low separability or because differentiation is unim¬ 
portant for the objectives of this study. Post-classification class merging only 
concerns non-vegetation classes. Separation of these classes has increased 
detection rate of bare ground and separation from artificial surfaces. Exam¬ 
ples of merged classes are: Alluvial fans. Unpaved ways and Mining area. 

Persistent shadows and clouds were removed from the classified maps using 
a customized focal filtering approach (modal-filter) with increasing window 
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Table 6: Description of the LULC classes used in this study. If major differences 
in class formulation exist additional descriptions are given. Classes 
merged to Bare Ground after classification are grey. 


La Campana 

CLASSES: 1-11, 20, 21 

(l)Mediterranean forest 

Natural forest vegetation. Evergreen, dense matorral, scle- 
rophyll forest, palm forest. Vegetation cover 75-100% 

(2)Deciduous forest 

Nothofagus forest (var. obliqua), occurrence only on El 
Roble, La Campana, Campanita summits (south/south-east 
facing slopes) 

(3)Shrubland 

Matorral Espinosa, 25-75% shrub cover ( Acacia chilensis 
etc.); evergreen; dominantly on north-facing slopes 

(4)Orchard 

Wine yards and fruit orchards (olives, avocado, oranges, 
mandarins); often on moderately steep slopes 

(5)Eucalyptus plantation 

Regular non-natural stands of Eucalyptus globulus and 
nitens with different age; evergreen, on moderate slopes 

(6)Pine plantation & 
Conifers 

Regular non-natural stands of Pirns radiata with different 
age; evergreen and coniferous, on moderate slopes 

(7)Grasslands 

Pasture, irrigated crops, “praderas” (i.e. ~ temperate grass¬ 
land) 

(8)Agriculture 

Bare agriculture fields (i.e. arable land), without growing 
crops at the date of image acquisition 

(9)Bare ground 

Uncovered ground (rocks, uncultivated soil and sediments) 

(lO)Alluvial fans 

Fans of rivers into lakes and before glaciers 

(1 l)Artificial surfaces 

Buildings, paved roads and other impervious human struc¬ 
tures (e.g. solar panels) 

Shadows (20), Water 
(21), Snow&Ice (22), 
Clouds (23) 


Nahuelbuta 

2, 5-8, 10-15, 21 


(2)Deciduous forest Natural forest (e.g. Nothofagus dombeyi), with less coni¬ 

fers/evergreen trees than deciduous species per pixel 


(12) Araucaria Araucaria araucaria, occurrence limited to Nahuelbuta 

National Park, especially on exposed summits, evergreen 

(13) Young growth Young growth natural vegetation and young growth planta¬ 

tions of both Eucalyptus and Pine 
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(14)Clearing 

Cleared areas (due to wood extraction, storms, forest 
fires),with nearly uncovered ground (i.e. with little under¬ 
story remaining) 

(15 )Unpaved ways 

Forest ways, higher brightness (devoid of vegetation) than 
“Clear cut”, elongated shape 

Santa Gracia 

9, 11, 16-23 

(16)Trees&tall shrubs 

Trees like Eucalyptus and Prosopis, no plantations. Ever¬ 
green, with high vegetation cover and height (>1.50 m). 
Situated in valleys, mostly close to roads or settlements 

(17)Shrabs&cacti 

Small shrub associations with cacti (e.g. members of the 
Lomas formation, Eulychnia, Nolana , Heliotropium, 
Tetragonia, or Euphorbia, etc.). Often in valleys or depres¬ 
sions and closer to the coast 

(18)Irrigated fields 

Crops grown due to watering, non natural grasslands 

(19)Mining area 

Remarkable mining areas (grey tone) (e.g. El Romeral) 

Pan De Azucar 

9, 10,17, 20, 21, 23 

Only vegetation class is “(17)Shrubs&Cacti” (e.g. 

Copiapoa cacti) 


size in R. For the purpose of improved visualization in non-scientific out¬ 
reach and calculation of object areas, rasters were converted to vector ge¬ 
ometry and maps were edited based on expert rales, e.g. all artificial surfaces 
inside the Ocoa catchment were reclassified to bare ground because no huge 
man-made structures are present. Further, all objects classified as orchard <= 
200 m 2 (i.e. 2 pixels) were reassigned to the class, which is closest judging 
from the J-M separability measures (chapter 4). As Araucaria are known to 
be limited in occurrence to the area inside Nahuelbuta National Park, all 
classified objects outside the park borders were merged to the general Coni¬ 
fer class. 

4.2.6 Machine Learning 

The classification procedure was done in R statistical programming software 
(R Development Core Team, 2008; Crawley, 2013) and two machine learn¬ 
ing approaches (Random Forest, Support Vector Machines) were compared 
regarding their suitability for LC classification. Programming in R has cer¬ 
tain advantages compared to traditional GUI based approaches: 
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1. The user is not bound to a certain software package 

2. Allows for automatization relevant for processing large amounts of data 

3. Improves transparency through better documentation and repeatability 
(Tippmann, 2015) 

A Random Forest consists of several decision trees aggregated to form a 
decision “forest”. This kind of ensemble learning can substantially improve 
predictive performance when compared to a classical simple decision tree 
(James et al., 2013). Random Forest is popular in remote sensing because it 
can handle high data dimensionality, multicollinearity and is robust against 
overfitting (Breiman, 2001; Belgiu and Dragut, 2016). In addition Random 
Forests often provide better land cover classification accuracies than for 
example the often used maximum likelihood approach (Waske and Braun, 
2009). Each tree in a Random Forest is built on a bootstrapped sample of the 
training data. At each split node of a single tree a random subset m of the 
predictors is evaluated. This randomness improves the simpler bagging by 
decorrelating single trees. The classification is reached by aggregation and 
majority voting of all decisions passed from individual trees (Liaw, A. and 
Wiener, M., 2002). With the simple core structure of a Random Forest, just 
two parameters (number of predictors at split m; number of trees n) need to 
be set. This easy parametrization is an advantage compared to other black¬ 
box machine learning approaches (Pelletier et al., 2016). 

As it is impossible to survey a whole ensemble of trees created by RF, varia¬ 
ble importance (VI) was assessed through the “varlmpPlot” package in R. In 
general, the variable importance plot shows, which variables form a major 
part of the predictive power of the model. If the uppermost variable is ig¬ 
nored and the model run a second time, the RMSE is expected to rise be¬ 
cause of the substantially lower predictive power. Same is the case for the 
other variables mentioned in the plot, while their contribution to the predic¬ 
tive power of the model is expected to be declining with the decrease in the 
position inside the plot (top down). Vice versa, variables unmentioned in the 
plot are potential candidates for omission in future versions of the model 
(James et al., 2013). As RF always takes bootstrapped samples to construct 
single trees, variable importance measures may slightly differ between dif¬ 
ferent mns of the model. Although Random Forests are said to be able to 
handle large numbers of predictor variables, data reduction based on variable 
importance has shown to increase classification accuracy (Belgiu and Dragu(, 
2016). However, no reduction of the feature space has been done in this 
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study, as initial classifications have not shown a gain in accuracy after omis¬ 
sion of low importance variables. 

SVMs are often used in remote sensing due to their ability to give good re¬ 
sults even with limited sample data, which is often an issue in remote sens¬ 
ing. According to their definition they are “supervised non-parametric statis¬ 
tical learning techniques” (Mountrakis et al., 2011). Various studies have 
been published describing the use of SVM classificators in all kinds of appli¬ 
cations from LULC to the mapping of different stages of imperviousness 
(Bachofer et al; Kavzoglu and Colkesen, 2009). Reviews are given by Huang 
et al. (2002), Mountrakis et al. (2011) and Pal and Mather (2005), with the 
latter showing a continuous increase of studies implementing SVMs in the 
field of remote sensing. SVMs work with so-called hyperplanes that try to 
linearly separate the data points (Vapnik, 1995). The borders between classes 
are chosen based on the size of the margin between the data points closest to 
the hyperplane (“Support vectors”). The “learning” takes place while finding 
the optimal boundary to separate the training samples. If linear separation is 
impossible, kernel functions help SVMs to project data into higher¬ 
dimensional feature spaces where the separating hyperplane can be con¬ 
structed. The configuration of the kernel and the kernel specific parameters is 
crucial but case dependent with little guidance by the literature (Pal and 
Mather, 2005). In this study, the chosen kernel is a radial basis function: 

K(xi,xj ) — exp(—g\\ xi - xj\\ 2 ),g > 0 (11) 

where x, and Xj represent the training data and class labels, and g is the gam¬ 
ma term in the kernel (Zakeri et al., 2017). Besides gamma, sigma and C (i.e. 
“Cost”) are important parameters of the SVM model. Sigma denotes the 
standard deviation of the underlying Gaussian distribution while the C-value 
defines the margin of the hyperplane and therewith how much the classifica¬ 
tion wants to avoid misclassification. Different sigma and Cost values were 
tested with gamma held constant using a gridsearch in R, which automatical¬ 
ly decides for the value resulting in the highest classification accuracy. 

4.2 .7 Accuracy Assessment 

Assessing the accuracy of produced maps is a central step in the classifica¬ 
tion process. It measures the quality on a statistical basis thereby increasing 
reliability. In this study, accuracy of every map was statistically quantified at 
least two times. Regarding the time series no accuracy measurement was 
conducted. The first step includes the cross-validation measures produced 
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during the model training process described in chapter 4.2.6. 10-fold cross- 
validation was chosen (Kohavi, 1995), which significantly decreased pro¬ 
cessing time compared to leave-one-out cross-validation (LOOCV). The 
model with the highest Kappa coefficient (Cohen, 1960) was finally taken for 
prediction. Cohen's Kappa coefficient was preferred to Overall Accuracy 
(OA) because it accommodates for the case that some pixels are allocated to 
the correct class purely by chance (Congalton and Green, 2009). Cohen's 
Kappa coefficient is calculated as follows: 


_ Po-Pc 
1 ~Pc 


( 12 ) 


with p 0 ='Z Pa being the sum of frequency in the diagonal of the error matrix 
and p c = Yj Pi + Pj+ is the hypothetical probability of chance agreement. 

In a second step, the independent validation datasets were taken to assess the 
accuracy of the produced map. Confusion matrices and a number of associat¬ 
ed measures such as producers, users and OA as well as the Kappa coeffi¬ 
cient were calculated. Coefficients >0.80 represent strong agreement and 
good accuracy, 0.40-0.80 is middle, <0.40 is poor agreement. A value of e.g. 
0.8 would represent an 80% higher agreement than expected by chance. Note 
that both OA, Kappa and especially class-wise statistics are presented in this 
study due to the recent criticism Kappa coefficient was exposed to regarding 
its suitability in measuring map accuracy through a single statistic (Pontius 
and Millones, 2011; Olofsson et al., 2014). Instead of users and producers 
accuracy (to which the respective R package refers as Precision and Recall), 
the FI-score is presented, which is a per class harmonic mean of both 
measures: 


„ „ 1 „ yrecision*recall .,,. 

F, = 2*—; -t-=2*- (13) 

-h- precision+recall 

recall precision 

The McNemar test (without correction for continuity) was applied to investi¬ 
gate, which pre-processing or machine learning algorithm performed signifi¬ 
cantly better in predicting land cover (Foody, 2004). McNemar's test is a 
nonparametric, always two-sided statistical test that can be applied to 2x2 
confusion matrices. In this case it was used to test if the differences in pro¬ 
portions of correctly and falsely classified samples between the two classifi¬ 
cations were statistically significant. 

Following the approach presented by Vanonckelen et al. (2013) the discrep¬ 
ancy in accuracy between solely atmospherically corrected image and the 
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image with the complete preprocessing chain applied was quantified by sam¬ 
pling the difference area of the two classifications produced based on these 
images. Hereby, the merit of topographic correction given by increase in 
accuracy should be quantified. This is done only for the coastal part of the La 
Campana region, which has the highest diversity in classes and strong topog¬ 
raphy (Coastal Cordillera). The procedure applied here can be outlined as 
follows: 
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After subtraction of both classifications and conversion to polygons, 50 
points were randomly generated inside each class defined by the classifica¬ 
tion based on topographically corrected data. Class affiliation was deter¬ 
mined by analysis of high-resolution imagery in GoogleEarth. The dataset is 
then split into different illumination zones based on a hillshade model of the 
study region. Finally, accuracy is determined separately for each illumination 
zone and classification. 

4.2.8 Time series Analysis 

Simple bi-temporal change detection methods such as image differencing 
and Change Vector Analysis (CVA) (Coppin et al., 2004), have limited ca¬ 
pabilities for detecting land cover changes within time series that are heavily 
influenced by seasonal climatic variations (Verbesselt et al., 2010). Long 
term trends are likely to be masked by seasonal variation and the results of 
the study are heavily dependent on the appropriate image selection (Watts 
and Laffan, 2014; DeVries et al., 2015). The BFAST (Breaks For Additive 
Seasonal and Trend) has proven successful in separating trend and seasonal 
components, while accurately detecting breaks in vegetation index time se¬ 
ries caused by natural and manmade disturbances (Verbesselt et al., 2010; 
Watts and Laffan, 2014; Schultz et al., 2016). BFAST is able to estimate 
time and number as well as magnitude and direction of changes (Verbesselt 
et al., 2010). (Verbesselt et al., 2010) have implemented their code in R sta¬ 
tistical programming software (R Development Core Team, 2008; Verbesselt 
et al., 2014). The code can be modified to tune the parameters “h” (minimal 
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segment size between potentially detected breaks in the trend mode) and 
number of iterations and maximum breaks to detect. 

Time range of the analysis was limited to 4 years in which Landsat 8 data 
were available (April 2013 till March 2017). This avoids the huge gap after 
the decommission of Landsat 5 in late 2011 and the launch of Landsat 8 in 
early 2013. Landsat 7 data after the malfunction in 2003 (SLC-off) were 
omitted due to the unavoidable striping, which cannot be fully removed from 
the images. NDVIs for BFAST analysis were calculated as means across the 
whole catchments. To fit the demands of the BFAST algorithm, missing data 
was imputed using spline interpolation (Moritz, 2017). Spline was preferred 
to linear interpolation because it was found to better reflect natural patterns. 
One must consider that interpolation may introduce artifacts in both frequen¬ 
cy and amplitude (Trauth et al., 2007). No smoothing was applied on the 
time series prior to BFAST analysis. Partly because significant outliers were 
already removed in advance. In addition, smoothing also has inherent disad¬ 
vantages by altering actually observed values (Pettorelli et al., 2005). 

To monitor changes on a more detailed level, height-specific and class- 
related time series were calculated as well. Therefore, the catchments were 
once separated into different altitudinal belts with 200m spacing and mean 
NDVIs calculated in every belt. Secondly, NDVIs are averaged across the 
areal coverage of each single class inside the catchments. The classification, 
which performed best was chosen to do so. Note that in the case of La Cam- 
pana, it was necessary to merge two independent classifications, because the 
catchment is located at the border of two scenes (s. chapter 4.1.1). No com¬ 
bined altitudinal and class differentiation was conducted. 
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This section first presents pre-classification results from the performance of a 
set of different topographic corrections and class separability measurements 
(chapters 4 and 4). Afterwards the outcomes of the classification using dif¬ 
ferent pre-processing and machine learning algorithms as well as the final 
maps are presented (chapters 5.3, 0 and 5.5). Chapter 5.6 deals with the three 
kinds of NDV1 time series calculated as part of this study, followed by a final 
summarizing interpretation of the most important results. 


5,1 Comparison of Topographic Corrections 

According to Table 7 the influence of topographic correction differs between 
bands. In general, the mean of the NIR band is subject to stronger alteration 
than the mean of bands from the visible spectrum. The sum of all means 
across the presented techniques is -0.053 in NIR, -0.01 in the Blue band and - 
0.001 in Green and Red. If also SD is considered, the Red band features the 
overall smallest deviation. The largest deviation from 0 in terms of mean 
reflectance value occurs after correction with the Sen2Cor algorithm. The 
negative sign implies that the mean after correction is higher than before, 
which is the case for most of the tested algorithms. A high deviation in mean 
does not necessarily involve a high deviation in SD. E.g. while deviations are 
high for both statistics if Cosine correction is used, the total change in SD 
resulting from Sen2Cor is relatively low. The top three algorithms, which are 
closest to 0 (considering both mean and SD) are as follows: 1. Minnaert, 2. 
C-correction and 3. C-NDVI. Only Minnaert shows perfect maintenance of 
both statistics. 

A subset of the presented correction results is depicted in Figure 8. During 
the time of image acquisition sun-caused shadowing is present on south¬ 
facing slopes, which cover the lower part of the images. Due to the different 
illumination, natural forests on south and north-facing slopes appear in dif¬ 
ferent colors of green. This color inhomogeneity is only partly explained by 
the occurrence of deciduous trees on south-facing slopes and thus differences 
in species composition and vegetation density. Visual comparison of the 
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Table 7: Comparison of 9 topographic correction techniques. Mean and Standard 
Deviation (SD) of the whole image are given for the 4 10m resolution 
bands. Total change is calculated with absolute values. All calculations 
based on the scene La Campana 1 (coastal image) and the TanDEM-X 
DEM. Abbreviations: “Impr. Cosine” = Improved Cosine. 




B 

G 

R 

NIR 

Total Change 

Cosine 

MEAN 

-0.008 

-0.007 

-0.005 

-0.013 

0.034 


SD 

-0.012 

-0.010 

-0.006 

-0.016 

0.045 

Impr. Cosine 

MEAN 

0.000 

0.000 

0.000 

-0.001 

0.001 


SD 

-0.003 

-0.002 

-0.001 

-0.006 

0.012 

scs 

MEAN 

-0.004 

-0.003 

-0.002 

-0.006 

0.015 


SD 

-0.007 

-0.006 

-0.003 

-0.010 

0.026 

Minnaert 

MEAN 

0.000 

0.000 

0.000 

0.000 

0.000 


SD 

0.000 

0.000 

0.000 

0.000 

0.000 

Minslope 

MEAN 

0.003 

0.003 

0.002 

0.006 

0.014 


SD 

0.000 

0.000 

0.000 

-0.001 

0.001 

C-correction 

MEAN 

0.001 

0.001 

-0.001 

0.001 

0.004 


SD 

0.001 

0.000 

0.000 

0.000 

0.002 

Gamma 

MEAN 

0.011 

0.010 

0.008 

0.018 

0.048 


SD 

0.001 

0.000 

0.001 

-0.003 

0.004 

C-NDVI 

MEAN 

-0.002 

-0.002 

-0.002 

-0.003 

0.009 


SD 

0.001 

0.001 

0.001 

0.001 

0.004 

Sen2Cor 

MEAN 

-0.011 

-0.001 

-0.001 

-0.055 

0.067 


SD 

0.000 

0.002 

0.001 

0.001 

0.004 


results from the three topographic corrections presented in Figure 8 shows 
that Sen2Cor performed best regarding color alignment. The pronounced 
green tone in the original image is adjusted to the brownish green of the well 
illuminated north-facing slopes and valley bottoms inside the Ocoa catch¬ 
ment. 
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71°5.940'W 71°5.<M0'W 71°5.940'W 71°5.040'W 


Figure 8: Comparison of topographic correction methods in the La Campana Na¬ 
tional Park. Clockwise from u.l.: Atmospherically corrected image, 
Sen2Cor, C-NDVI, Minnaert. The arrow marks the approximate position 
at the time the image was taken (280° solar azimuth at 14:37 o'clock lo¬ 
cal time). 


C-NDVI and Minnaert are closer to the original image preserving the imbal¬ 
anced colors. Unvegetated spots on both north and south-facing slopes are 
occasionally overcorrected resulting in very bright colors in the Sen2Cor and 
Minnaert algorithms, while C-NDVI seems to attenuate overall differences in 
brightness. All algorithms except from Gamma and Improved Cosine have in 
common that in areas of strong topography with shadowing (e.g. scarps on 
Cerro La Campana, not shown in Figure 8), artifacts in strange colors (blue- 
green) occur. 
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5.2 Class Separability 

Separability of classes was tested using Jefffies-Matusita (J-M) distance. 
Values >1 show tendency for better separability, while classes with <1 are 
potentially less separable thus increasing the likeliness of classification er¬ 
rors. The results of the assessment for Pan De Azucar are shown in Table 8. 
For other study areas, only exemplary values are presented and the reader is 
referred to the supplementary material for a full list. 

For Pan De Azucar in general, most classes have at least one class from 
which separability is below the mean of 1. Only water is well distinguished 
from every other class, with a maximum separability of 1.35. The overall 
mean across all 6 classes in Pan De Azucar is 0.94. The lowest mean separa¬ 
bility features the bare ground class (0.71), while water has the highest 
(1.25). Low separability for vegetation (Shrubs&Cacti) exists with the bare 
ground class (0.37), which is the lowest separability measure of all classes in 
Pan De Azucar. Apart from the case mentioned before separability does not 
fall below a critical value of 0.5. 

Between the 10 initial (before merging of artificial surfaces and mining areas 
with bare ground) classes in Santa Gracia overall separability is lower than in 
the other study areas (0.84). Mean class separability both among the 3 vege- 


Table 8: Separability of classes for training samples acquired in Pan De Azucar. 

Classes with potentially low separability and values smaller than a 
threshold of 0.8 are highlighted in grey. Abbreviations: “ArtSurf ’ = Arti- 
cial Surfaces. 



ArtSurf 

Bare 

Cloud 

Shadow 

Shrub- 

Cacti 

Water 

ArtSurf 

0 

0.57 

1.03 

1.02 

0.78 

1.12 

Bare 

0.57 

0 

0.77 

0.66 

0.37 

1.18 

Cloud 

1.03 

0.77 

0 

0.86 

1.03 

1.35 

Shadow 

1.02 

0.66 

0.86 

0 

0.81 

1.29 

ShrubCacti 

0.78 

0.37 

1.03 

0.81 

0 

1.31 

Water 

1.12 

1.18 

1.35 

1.29 

1.31 

0 

Mean 

0.90 

0.71 

1.01 

0.93 

0.86 

1.25 
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tated classes and all classes is lowest for shrubs&cacti (0.71). Better differen¬ 
tiation exists for trees&shrubs (0.87) and irrigated fields (0.95). Only 
snow&ice exceeds a mean value of 1 (1.11). Regarding intra-class separabil¬ 
ity, the lowest values exist between shrubs&cacti and trees&shrubs (0.31), 
while J-M between shmbs&cacti and irrigated fields is 0.7 or 0.58 between 
irrigated fields and trees&shrubs. 

Overall mean J-M in La Campana 1 is 0.89 and 1.1 in La Campana 2. Mean 
vegetation class separabilities are generally higher in La Campana 2(1: 0.93; 
2: 1.08; 3: 1.00; 4: 1.14; 5: 1.20) than in La Campana 1 (1: 0.74; 2: 0.94; 3: 
0.88; 4: 0.87; 5: 0.69; 6: 0.85). Eucalyptus plantations have the lowest values 
amongst all classes in La Campana 1 (0.69), but not in La Campana 2 (1.20). 
Separability from Eucalyptus is low for natural forests and pine plantations 
(1: 0.31; 6: 0.35). Plantations of pine are missing in La Campana 2. 

In Nahuelbuta, none of the vegetation classes (2: 0.72; 5: 0.81; 6: 0.71; 7: 
0.77; 12: 0.82; 13: 0.82; 14: 0.91) is above the study area mean (0.94). While 
the J-M value between plantation types (classes 5&6) is only 0.29, separabil¬ 
ity of clear cut from other vegetation classes is always above 0.85. Similar to 
the other study areas water surfaces show generally high measures of separa¬ 
bility (1.40) and artificial surfaces as well as agricultural fields separate less 
well from each other (0.47). 


5.3 Classification Accuracy 

Two kinds of machine learning algorithms were used in this study and their 
performance was tested inside each study area. Results of this assessment are 
presented in chapter 5.3.1. In La Campana 1 a RF classification was per¬ 
formed three times, once with the data corrected for topography using the 
SRTM DEM, a second time with the TanDEM-X DEM and finally without 
topographic correction (chapter 5.3.2). 

5.3.1 Impact of the Classifier 

Table 9 shows that RF outperformed SVM in every study area, although the 
differences in overall classification accuracies are very small. The hugest 
gain in accuracy through RF exists in La Campana 2, where the difference in 
OA amounts to 0.04, which is still just a minor increase. Significance testing 
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Table 9: Comparison between RF and SVM in La Campana 1 (OA/Kappa Coeffi¬ 
cient). Significance of the difference tested with McNemar's test 
(McNemar's Chi-squared/p-value) and 1 degree of freedom. 


Study Area 


RF 

SVM 


McNemar 

Pan De Azucar 


0.99/0.99 

0.97/0.96 


2.41/0.12 

Santa Gracia 


0.98/0.98 

0.96/0.96 


1.29/0.26 

La Campana 1 


0.87/0.85 

0.85/0.84 


1.91/0.17 

La Campana 2 


0.86/0.84 

0.82/0.79 


6.49/0.01 

Nahuelbuta 


0.90/0.89 

0.88/0.86 


2.40/0.12 

with McNemar's test on the 5% level shows that only in La Campana 2 the 
gain in accuracy with RF is statistically significant. Complete omission of 
the Agriculture class with SVM is the major reason for this discrepancy 
(Table 10). However, SVM did not perform better than RF regarding the 
other classes except for Shadows. 

Table 10: Class-wise FI 

-score for RF and SVM in La Campana 2. 


Class 

RF 

SVM 

Class 

RF 

SVM 

1 

0.79 

0.77 

8 

0.53 

0 

2 

1.00 

0.98 

9 

0.81 

0.77 

3 

0.72 

0.68 

11 

0.85 

0.90 

4 

0.90 

0.83 

20 

0.83 

0.90 

5 

1.00 

0.99 

21 

0.90 

0.89 

7 

0.96 

0.93 

22 

1.00 

0.99 


Although OA in La Campana 1 suggests a slightly better performance of RF, 
7 classes were predicted better with SVM and only 5 with RF. Albeit, in 9 
cases differences are smaller than 0.02. Again, SVM did better in delineating 
shadows (RF: 0.36; SVM: 0.48) but had difficulties with shrubland, which 
was often confused with natural forest or bare ground (RF: 0.92; SVM: 
0.76). In Nahuelbuta, RF outperformed SVM in 7 of 11 classes. The highest 
plus was in favor of SVM (Araucaria, RF: 0.95; SVM: 0.99), which is close¬ 
ly followed by RF and the Eucalyptus class (RF: 0.84; SVM: 0.81). 
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In Pan De Azucar and Santa Gracia, variation between classifiers is negligi¬ 
ble. While in Santa Gracia RF did better in every case with the maximum 
accuracy plus being 0.03 (Bare), SVM was slightly better in detecting vege¬ 
tation in Pan De Azucar (0.05) but lagged behind in every other case. 

5.3.2 Impact of Topographic Correction and DEM 

The left part of Figure 9 shows a visual comparison of TanDEM-X and 
SRTM DEM on the western part of the La Campana National Park including 
the Cerro La Campana summit. The degree of detail within the TanDEM-X 
DEM is higher, which is most obvious considering the ridges and uneven 
terrain radiating from the Cerro La Campana. The right part of the figure 
shows the remainder if each two DEMs are subtracted from another. Three 
different DEMs were used for the subtraction, the TanDEM-X DEM with 12 
m resolution, a 30 m SRTM DEM and a Lidar DEM with 2 m resolution, 
which has been generated within the EarthShape project. The elevation val¬ 
ues of the Lidar DEM are assumed to be most accurate and closest to true 
height a.s.l. and therefore serve as reference. For mathematical operations, 
DEMs were resampled to 10 m resolution. The range between Min and Max 
is higher if SRTM was subtracted from the Lidar DEM and lower for the 
TanDEM-X DEM, while the mean difference is huger for the TanDEM-X 
DEM (31.39 m; SRTM - Lidar: 10.40 m). The mean difference between 
TanDEM-X DEM and SRTM amounts to 21.39 m. In contrast to the SRTM, 
the difference of the TanDEM-X DEM resembles the terrain with the lowest 
offset on east, south-east facing slopes. The difference in the SRTM DEM is 
lowest on ridges but aspect dependency is less pronounced. Both TanDEM-X 
and SRTM DEM feature the highest difference on the bottom of incised 
valleys and at the steep eastern mountainside of the Cerro La Campana. 

Statistical figures do not show a significant difference in classification accu¬ 
racy if the TanDEM-X DEM or the SRTM DEM is used for topographically 
correcting both optical and radar data (Table 11). Class-wise comparison 
reveals that the TanDEM-X DEM data processing leads to an improved pre¬ 
diction of 3 vegetation classes including classes 1, 2 and 7. In case of the 
other 

vegetation as well as most of the remaining non-vegetated classes, pro¬ 
cessing with SRTM DEM leads to a lower error rate, whereby most of the 
differences are negligible and below 0.01 in FI-score. Although still statisti- 



48 


5 Results 



Figure 9: Comparison of the 30 m-SRTM (mid picture) and 12 m-TanDEM-X 
DEM (lower picture) with a Lidar DEM in the La Campana National 
Park. Small picture shows the position of the scene inside the park. The 
details visible in the TanDEM-X are higher compared to the SRTM espe¬ 
cially in case of the small-scale topography surrounding the Cerro La 
Campana (summit on the left). Images on the right show the differences 
if two DEMs are subtracted from another. Color scale on both DEMs and 
differences range from black (low values) to white (high values). 


cally insignificant, the pure atmospheric correction (in the following referred 
to as “Atmo”) has the highest OA statistics. Similar to the processing with 
SRTM DEM, the TanDEM-X DEM did better in detecting classes 1, 2 and 7, 
whereby the differences are smaller compared to SRTM except for class 7. In 
case of vegetation classes 4, 5 and 7 the SRTM DEM was slightly above the 
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Table 11: Comparison between classifications based upon data corrected with 
TanDEM-X, SRTM and without DEM (OA/Kappa Coefficient). Signifi¬ 
cance of the difference tested with McNemar's test (McNemar's Chi- 
squared/p-value) always against the TanDEM-X DEM processing and 1 
degree of freedom. 


Study Area 

La Campana 1 

McNemar 

TanDEM-X 

0.87/0.85 

- 

SRTM 

0.87/0.86 

0.03/0.86 

NONE 

0.90/0.89 

0.28/0.60 


atmospheric correction. The main accuracy gain achieved by atmospheric 
correction was within the Shadow class, which does not form part of this 
study's objectives. After exclusion of Shadows from the calculation of the 
mean FI-score, the results are still favoring atmospheric correction, but the 
difference to correction with TanDEM-X is just 0.01 (0.03 with the Shadow 
class). 

Figure 11 shows the different illumination zones derived from the TanDEM- 
X DEM under the lighting conditions on acquisition time. North and north¬ 
western exposed slopes are well illuminated while southerly facing slopes 
receive little light. A total of 600 random points was generated in the differ¬ 
ence area between both classifications sampling all three illumination zones. 
All points belonging to the Shadow class were removed by what the final 
number of points was 546. The affiliation to illumination zones among sam¬ 
ple points was as follows: 77 in low, 220 in moderate and 249 in the highly 
illuminated zone. Both classifications disagreed in a total area of 150.4 km 2 . 
Regarding areal coverage, the top three in mismatched classes were: 

• 3 (topographically corrected, in the following referred to as “Topo”) -> 7 
(Atmo) = 15.23 km 2 

• 7 (Topo) 4 (Atmo) = 11.69 km 2 

• 11 (Topo) 8 (Atmo) = 10.68 km 2 

In the low and high illumination zones, atmospherically corrected data leads 
to better results in OA than topographically corrected (Figure 10). Fligher Fl- 
scores are recorded for classes 1 (Topo: 0; Atmo: 0.51) and 3 (0.26; 0.31), 
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0.40 

0.30 

0.20 

0.10 

0.00 


Figure 10: OA in the difference area between RF classifications based on atmos¬ 
pherically and topographically corrected data. “Topo” = atmospherically 
and topographically corrected data; “Atmo” = “atmospherically corrected 
data”. 

but especially in class 1 the discrepancy is remarkable because in Topo no 
Natural forests at all agree with the reference. Each 6 Natural forest pixels of 
20 in total were wrongly assigned to classes 2 and 4 and 4 pixels to class 5, 
which has the lowest separability according to J-M. Classification with At¬ 
mo assigned 16 out of 20 Natural forest pixels correctly, whereby the still 
low FI-score is explained by a low Producers Accuracy (0.37). Classes 3 
(0.57; 0.25) and 9 (0.28; 0.06) are predicted with higher accuracy within 



Low Topo Low Atmo Moderate Moderate High Topo High Atmo 
Topo Atmo 



Figure 11: Illumination zones in La Campana with a solar azimuth of 73.1. 
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Topo. Within the moderately illuminated zone, Topo performed better espe¬ 
cially in case of Water (0.44; 0) and Pine (0.38; 0). Water (0.50; 0.13) was 
also more accurately predicted in the high illumination zone, while detection 
of Shrubland (0.26; 0.54) and Eucalyptus (0; 0.29) was better within Atmo. 
All in all, 15 times out of 33 Atmo achieved higher FI-scores than Topo, in 
comparison to 10 times in which Topo did better than Atmo. Albeit, if slight 
discrepancies are excluded from the data and only differences >0.2 are con¬ 
sidered the relation is in support of Topo (8:3). 


5.4 Final LULC Maps 

Based upon the results presented in the previous chapters RF classifications 
were chosen for the final maps and SVM classifications were rejected. Huger 
maps can be found in the supplementary material. LULC frequency tables 
presented in this chapter represent the original state before reclassification 
and after shadow removal. LULC maps presented herein are reclassified as 
specified in the text. The modal filter described before was applied to Pan De 
Azucar, Santa Gracia and La Campana. Only in Pan De Azucar and Santa 
Gracia also Clouds were present (Figure 12), whereby Shadows in the other 
study regions were solely topographic Shadows. Clouds and Shadows in Pan 
De Azucar cover together an area of 10.76 km 2 . The final maps show a mis- 
classification of Artificial surfaces north of the catchment and confusion with 
vegetation along a sandy road (C-l 10), which first runs from the coast to the 
catchment where it splits into the C-l 10, which crosses through the catch¬ 
ment and the C-l 12, which continues in northwestern direction. Vegetation 
is concentrated on the coast and practically none was classified further than 
50 km from the coast. This is supported by the distribution of classes with 
height (Figure 13). Thus, plants cover less than 1% inside the whole study 
area. Within the catchment percent cover is about 2%, which is remarkably 
enough higher than the relative amount covered by plants inside the Santa 
Gracia catchment. The hugest connected patch of class 16 is located close to 
the coast and measures 18 000 m 2 , while the majority of vegetation objects is 
very small only 2% exceed 1000 m 2 in size. Two highways cross the Pan De 
Azucar region, one about 20 km from the coast and one further inlands (~65 
km). Except from the parts obscured by clouds, the coastal highway is clearly 
visible. Water is scarce in both Pan De Azucar (Table 12) and Santa Gracia 
(Table 13), while in the latter the only noteworthy inland water is a salt lake 
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Figure 12: RF classification map for the Pan De Azucar region. The upper map 
shows the classification before Shadow&Cloud removal. In the lower 
map, Shadows&Clouds are removed and the outlines of the central 
catchment highlighted in black and background in the same color used 
for Water. 


Table 12: Distribution of classes inside the whole Pan De Azucar region (in km 2 
and as % of total area) and in the catchment (“Area C.”). Shadows and 
clouds are removed from the data. 


Class 

Area [km 2 ] 

% 

Area C. [km 2 ] 

% 

9 

1293.08 

98.94 

22.31 

97.44 

11 

2.37 

0.18 

0.04 

0.17 

17 

11.23 

0.86 

0.55 

2.39 

21 

0.25 

0.02 
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Height a.s.l. [m] 


Figure 13: Cumulative distribution of LULC in the Pan De Azucar region across 
different height zones (200 m steps) in % of total area inside each height 
zone. The lowest height zone always includes heights between 0 and 200. 
Due to the dominance of class 9 the y-axis focusses on values between 
90-100%. 


Table 13: Distribution of classes inside the whole Santa Gracia region (in km 2 and 
as % of total area) and in the catchment ("Area C."). Shadows are re¬ 
moved from the data. 


Class 

Area [km 2 ] 

% 

Area C. [km 2 ] 

% 

9 

2395.95 

94.34 

80.92 

98.67 

16 

1.43 

0.06 

0.01 

0.01 

17 

134.52 

5.30 

1.08 

1.32 

18 

0.85 

0.03 



21 

0.19 

0.01 



22 

6.82 

0.27 




~55km from the coast (0.23km 2 ), in Santa Gracia the vector map records 194 
independent water bodies. These are partially only few pixels in size and the 
hugest is roughly 60 000 km 2 . Percent cover of water surfaces in Pan De 
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Azucar is so small that it is barely visible in Figure 13. The maximum con¬ 
tribution to LULC is between 800-1000 m (0.04%), at which the salt lake 
occurs. 

Three plant communities are distinguished in Santa Gracia. Regarding areal 
extension in both the whole study area as well as inside the catchment, class 
17 is at least 94 times more frequent than classes 16 and 18. 

Because of its dependency on irrigation, class 18 is missing inside the 
catchment. As mentioned before, the relative dominance of bare ground is 
even greater than within the Pan De Azucar catchment closer to the desert 
core. Higher vegetation from class 16 inside the catchment covers 0.01 km 2 , 
which equals 100 pixels. In the whole region, patches belonging to class 16 
are generally small, the hugest covers 0.02 km 2 . 


7ri2.0(yw 71”0.00'W ZOHS.OO'W 70'36.00'W 70"24.00'W 70U2.00 W 



Figure 14: RF classification map for the Santa Gracia region after Shadow removal. 

Catchment outlined in black. 

Close to the coast, there is a misclassification within shadowed south ex¬ 
posed slopes, which are spectrally close to vegetation pixels. Irrigated fields 
exist close to the coast north of La Serena, which itself is outside the study 
area, and further inlands in secondary valleys of the main Elqui valley. The 
hugest irrigated field classified as one object measures 211 900 m 2 , the larg¬ 
est overall vegetation object belongs to class 17 and covers 16.25 km 2 in an 
assemblage of loosely connected Shrubs&Cacti located close to the coast. 
Apart from the major occurrence close to the coast, class 17 is detected ~60 
km from the coast along an up to 3800 m high mountain ridge (cp. Figure 
15). In some cases, shadows at this place are again misclassified as higher 
vegetation. Topographic shadowing increases with height so that at areas 
>3000 m 1.64% information is lost due to shadowing, while shadows thrown 
by clouds account to maximum 0.15% (<200 m). Snow & ice in the summer 
image first occur at heights >2800 m. 
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Figure 15: Cumulative distribution of LULC in the Santa Gracia region across dif¬ 
ferent height zones (200 m steps) in % of total area inside each height 
zone. Due to the dominance of class 9 the y-axis focusses on values be¬ 
tween 80-100%. 

Extent of classes in the whole La Campana region and inside the Ocoa 
catchment is presented in Table 14. Orchards that fall below a minimum 
object size of 300 m 2 were subsequently merged with the Natural forest class, 
to which spectral separability was lowest in both La Campana 1 and La 
Campana 2. Altogether, natural vegetation (classes 1, 2 and 3) other than 
Grassland, where a differentiation between natural grassland and pasture was 
impossible, amounts to 1114 km 2 . Plantations and orchards, which together 
form the non-natural vegetation, account for 82.3 km 2 , whereto Plantations 
contribute with 59.2 km 2 . Plantations, agricultural fields as well as settle¬ 
ments are more frequent in lowlands <600 m (Figure 18 & Figure 19) and 
reduce in frequency at higher altitudes. Both in La Campana 1&2, the rela¬ 
tive abundancy of shrubs and Mediterranean forest decreases rapidly within 
altitudes higher than 1400-1600 m. The relative dominance of bare ground 
especially in highly elevated areas is due to the low vegetation coverage in 
the high mountains of the Andes. Except from vegetation close to the Rio 
Aconcagua and its main tributaries, the Andes are largely devoid of vegeta¬ 
tion, which could be detected within 10 m resolution. 
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Table 14: Distribution of classes inside the whole La Campana region (in km 2 and 
as % of total area) and in the Ocoa catchment ("Area C."). Shadows are 
removed from the data. 


Class 

Area [km 2 ] 

% 

Area C. [km 2 ] 

% 

I 

416.04 

14.86 

16.73 

34.16 

2 

27.52 

0.98 

3.91 

7.99 

3 

670.60 

23.96 

21.77 

44.46 

4 

23.11 

0.83 

0.03 

0.05 

5 

50.87 

1.82 

0.04 

0.08 

6 

8.31 

0.30 

0.01 

0.02 

7 

152.81 

5.46 

0.01 

0.02 

8 

110.31 

3.94 

0.06 

0.12 

9 

1189.51 

42.49 

6.22 

12.69 

11 

95.69 

3.42 

0.20 

0.41 

21 

10.11 

0.36 



22 

44.47 

1.59 




In the whole study region out of the vegetation classes, Shrabland has the 
largest areal extent, while Orchards are spatially restricted to areas 
north/northwest and small occurrences east/northeast of the catchment. The 
largest connected communities of class 3 have been detected north of the 
Huechun reservoir (63.5 km 2 ) in La Campana 2 and north of the Ocoa 
catchment (87.1 km 2 ) in La Campana 1, with the latter being huge but very 
loosely connected. This one patch in La Campana 1 constitutes alone for 
~25% of all patches belonging to class 3 in this area. Whereas the contribu¬ 
tion of the largest Shrabland patch in La Campana 2 is 19.4%. Major planta¬ 
tions of Pine have been classified close to the coast and in the Valley sur¬ 
rounding Limache. Eucalyptus plantations are spatially 6 times more fre¬ 
quent than Pine plantations, whereby no Pine plantations were present in La 
Campana 2. Deciduous trees dominate in high elevation areas of the southern 
ridge connecting the peaks of El Roble and Cerro La Campana. Most class 2 
patches are rather small, the hugest forest measures 3.6 km 2 , while 23 894 
out of 24 350 patches do not fulfill the FAO standards of a forest (>0.5 ha). 
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3 - Shrubland 
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I I Ocoa Catchment outline 


Figure 16: RF classification map of La Campana 1 after Shadow removal. Ocoa 
catchment outlined in black and background in the same color used for 
Water (Pacific Ocean). 


Note that patches of single classes are considered individually and that the 
amount of forest patches will be potentially higher if all forest classes (and 
thus patches) will be merged. 



LULC classes 

■ 1 - Mediterranean Forest [ _ 5 - Eucalyptus Hi 11 - Artificial Surfaces Background 

I_ 2 - Deciduous Forest I 7 - Grass/Pasture B I 21 - Water I I Catchment outline 
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■ 4 - Orchards 9 - Bare Ground 0 


Figure 17: RF classification map of La Campana 2 after Shadow removal. Ocoa 
catchment outlined in black. 
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Figure 18: Cumulative distribution of LULC in the La Campana 1 region across 
different height zones (200 m steps) in % of total area inside each height 
zone. 



Figure 19: Cumulative distribution of LULC in the La Campana 2 region across 
different height zones (200 m steps) in % of total area inside each height 


zone. 
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Inside the Ocoa catchment, Mediterranean Forest and Shrubland dominate in 
terms of vegetation. Non-natural vegetation accounts for only 0.15% of the 
total surface area. The fact that Orchards are found in the catchment is due to 
the catchment borders, which differ from the borders of the national park. 
Therefore, to the north some fruit plantations are included. While at some 
places, e.g. northwest of the catchment, there is almost no visible edge be¬ 
tween both scene classifications, inside the catchment itself the difference is 
more obvious. This is due to the fact that areas sparsely covered by shrubs in 
highly elevated areas close to El Roble are classified as Shrubland in La 
Campana 1 but classified as bare ground in La Campana 2. Furthermore, the 
classifier in La Campana 2 tended to overclassify Mediterranean forests in¬ 
side the catchment, while at the same place Deciduous forests are overclassi¬ 
fied in La Campana 1. Deciduous forests close to El Roble on the southern 
ridge of the catchment are perfectly aligned between both classifications, 
while other land cover mapped at this place shows sharp edges (classes 1, 9 
and 11). 14.2% of all Deciduous trees fall inside the catchments boundaries. 
If all features that are close to or intersect the catchment are included into 
this calculation, 44.8% of all deciduous trees are addressed. 

Figure 20 shows the highly fragmented landscape in the Nahuelbuta region. 
Agriculture is frequent east and west of the Cordillera de Nahuelbuta. Close 
to the coast extensive Pine plantations exist (Max. = 287 600 m 2 ), while the 
westernmost part of the study region is dominated by large Eucalyptus plan¬ 
tations (Max. = 4.32 km 2 ). The cordillera itself is mostly forested with a 
mixture of natural growth and planted forests as well as wood extraction 
sites. In total, plants (classes 2, 5, 6, 7, 12, 13) cover 64.85% inside the study 
region (Table 15). Deciduous forests are by far the most frequent of all vege¬ 
tation classes covering more than 14 of the whole study region and more than 
Vi of the surface area inside the catchment. Consequently, the largest vegeta¬ 
tion patch belongs to class 2 and is located inside the Coastal Cordillera 
measuring 143.39 km 2 . 

The classification mapped 14.026 vegetation patches belonging to classes 2, 
5, 6 or 12 that fulfill the FAO criterion of a forest (>0.5 ha). Thereof 673.15 
km 2 belong to class 2 (90.58% of the total in the study region), 181.35 km 2 to 
class 5 (86.37%), 216.63 km 2 to class 6 (99,39%) and 17.13 km 2 to class 12 
(22.20%). This shows that while Eucalyptus is often found in large planta¬ 
tions, connected Araucaria forests are rare. On the coastal plains classes 
8 and 11 have been wrongly assigned to agriculture and settlements are 
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5 Results 


Table 15: Distribution of classes inside the Nahuelbuta region (in km 2 and as % of 
total area) and in the catchment ("Area C."). 


Class 

Area [km 2 ] 

% 

Area C. [km 2 ] 

% 

2 

743.18 

26.64 

15.82 

55.18 

5 

209.97 

7.53 

0.21 

0.73 

6 

217.96 

7.81 

0.15 

0.52 

7 

228.17 

8.18 

0.32 

1.13 

8 

487.36 

17.47 

0.07 

0.23 

9 

130.80 

4.69 

0.56 

1.96 

11 

54.54 

1.95 

0.10 

0.36 

12 

77.15 

2.76 

5.71 

19.92 

13 

332.95 

11.93 

2.47 

8.61 

14 

301.12 

10.79 

3.26 

11.37 

21 

6.99 

0.25 




missing as actual land use. Patches belonging to those classes situated on the 
coastal planes are therefore reclassified to the spectrally closest class (class 9, 
bare ground). Lakes west of the Cordillera are characteristically fringed by 
Conifers. River incised valleys that show a predominant northeast to south¬ 
west direction are, in contrast to the surrounding countryside, not used for 
agriculture. Extensive Grass/Pasture as well as forests of mainly deciduous 
species cover those areas. In general, the eastern area and other lowlands are 
dominated by settlements and associated land use (agricultural fields, pas¬ 
ture; Figure 21). Their abundancy decreases with height, while deciduous 
forests constantly increase until Araucaria and open mostly natural growth 
types (classes 13 & 14) take over their place. Riparian vegetation in the 
westernmost part of the study area consists of deciduous species forming 
narrow forest stands that rarely exceed 150 m width. All of these rivers are 
nearly completely covered with vegetation so that patches of water surface 
are seldom visible. Two huge patches of circular irrigation assigned to class 
7 are present in the study region west and southwest of Mulchen. The smaller 
of the two, which is located further south measures 679 500 m 2 , the other 900 
600 m 2 with a diameter close to 1 km. Misclassification in this part of the 
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Figure 20: RF classification map of the Nahuelbuta region. Catchment outlined in 
black. 

study region exists between agricultural fields and artificial surfaces. Class 
11 patches clearly identified as wrongly classified are assigned to class 8 
instead. 

Apart from class 2, which is by far the most dominant vegetation class, clas¬ 
ses 12 and 13 together account for 28.53% of the total vegetation cover in¬ 
side the catchment. Plants cover 86.09% of the surface area in the catchment 
if the mixture class 14 is ignored. Otherwise 97.46% are covered. Although 
only 7% of the overall Araucaria patch area is classified inside the catch¬ 
ment, 23% intersect the catchment including the largest Araucaria patch in 
the whole study region (12.38 km 2 ). More than 85% of all class 12 pixels are 
located on elevations higher than 1000 m and 80% on slopes steeper than 
10°. All patches that did not intersect Nahuelbuta National Park are reclassi¬ 
fied to class 6. Huge open areas (class 14) exist on less inclined surfaces 
(<1%) with more than 1200 m elevation in the eastern part of the catchment. 
In general, within the eastern part open areas and young growth are more 
frequent than within the central catchment. Here vegetation cover is dense 
and extensive forests of deciduous trees and smaller patches of Araucaria 
exist. Open areas again increasingly occur to the west (classes 7 & 14). The 
percentage of explicitly anthropogenic and non-natural classes inside the 
catchment is low (0.59%; classes 8 and 11). Artificial surfaces in the eastern 
part of the catchment are due to confusion with bright bare ground pixels on 
steep slopes and are reclassified to class 9. Agricultural fields in the catch¬ 
ment are reclassified to class 14, which is closest according to separability 


measures. 
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Figure 21: Cumulative distribution of LULC in the Nahuelbuta region across differ¬ 
ent height zones (200 m steps) in % of total area inside each height zone. 


5.5 Variable Importance 

This chapter presents the results of VI measures. A vegetation mask was 
applied to the study sites with high plant cover (La Campana, Nahuelbuta), 
so that VI can be given based on all or vegetation classes only. Furthermore, 
VI can always be given class-wise emphasizing the significance of single 
variables for the classification of certain classes. 

Great differences exist between the variables considered most important 
between the study regions, however all top 15 rankings have in common that 
radar variables are missing (Figure 22 & 23. Three variables are present in all 
top rankings considering all variables, these are RedEdge 3, NIR and 
GNDVI, but only the 3 rd Red Edge and the NIR band are likewise always 
inside the top 10. In general, bands from the near infrared region are fre¬ 
quently associated with high MDA values. Figure 22 and 23 show 20 of 24 
possible cases that each one of the 4 Red Edge or the NIR band could occur 
in the top 15. Of these 5 bands, the 2 nd Red Edge band has the lowest im¬ 
portance score in case a mean is calculated over the study areas (3.41 %). 
The five highest ranked variables are: 1. DEM (5.56%), 2. Red Edge 4 
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Variable Importance: Pan De Azucar 
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Variable Importance: Santa Gracia 
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Variable Importance: La Campana 2 
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Figure 22: VI for the whole set of classes listed by MDA. 

(4.76%), 3. Red Edge 3 (4.43%), 4. NIR (4.34%) and 5. SWIR I (4.18%). 
The five least ranked variables are: 1. Plan Curvature (0.05%) 2. Profile Cur¬ 
vature (0.07%) 3. Aspect (0.12%) 4. Maximum Probability (0.14%) 5. 
GLCM Correlation (0.17%). Of the four groups of variables (optical bands, 
indices, radar variables and DEM features) from all five top 15, 34 are indi¬ 
ces, 32 are optical bands and 8 are DEM features. Every index and every 
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Variable Importance: Nahuelbuta 


SWIRl 

5.82 

RedEdge2 

5.66 1 

SWIR2 

4.93 | 

DEM 

3.73| 

SLA VI 

3.66 | 

NIR 

3.51 | 

RedEdge3 

3.44 | 

Green 

3.41 I 

RedEdge4 

3.27 | 

CTVI 

3.23 | 

MSAVI2 

3.13| 

SWI 

3.11 | 

Blue 

3.07 | 

AutNdvi 

3.05 1 

GNDVI 

3.05 | 


0 1 2 3 4 5 6 

Mean Decrease in Accuracy [%] 

Figure 23: VI for the whole set of classes in Nahuelbuta listed by MDA. 

optical band except for the Green band appear in the top 15, while 8 DEM 
features are absent. Regarding the individual study sites, elevation is present 
in every ranking except for Pan De Azucar. Here, only optical bands and 
indices derived from optical bands can be found. DEM and derived features 
are considered most important in La Campana 2 where DEM and Positive 
Openness sum up to 12.06%, followed by Santa Gracia where DEM and TRI 
sum up to 11.11% of the total predictive power. The remainder is approxi¬ 
mately equally shared between optical bands and derived indices. Importance 
of seasonality information expressed by the seasonal NDVI (Spring, Au¬ 
tumn, Winter) increases from Pan De Azucar (Sum: 4.23%) to La Campana 2 
(17.10%) and decreases afterwards (Nahuelbuta sum: 8.14%). 

Figure 24 depicts the VI for the RF classification masking out all non¬ 
vegetation classes. The order of variables in the rankings differs from that 
derived from all classes. Between both VI estimations the hugest discrepan¬ 
cies in single importance scores exist for DEM, IRECI and GEMI, whereby 
all of them show higher importance with VI estimation using vegetation 
classes only. The relative importance of the DEM in La Campana 1 and in 
general is even higher than for VI with all classes (Mean: 8.80%). Just in 
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Variable Importance: La Campana 1 

Variable Importance: La Campana 2 
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Figure 24: VI for the vegetation classes only and the Araucaria class (l.r.). 

Nahuelbuta the DEM does not have the highest VI but is still on the 2 nd 
place. Though, DEM derived features other than DEM are present only in La 
Campana 1, so that the overall contribution of this feature group to the im¬ 
portance scores is lower than that of the infrared bands. The five highest 
ranked variables are: 1. DEM, 2. RedEdge4 (5.87%), 3. NIR (5.18%), 4. 
GEMI (4.19%), 5. RedEdge3 (4.15%). Overall the composition of variables 
presented in Figure 24 is: 21 optical bands, 19 spectral indices, 5 DEM de- 
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rived features. The newly established S-2 vegetation indices (IRECI, S2REP, 
verrelst) are all inside the 25 most important predictor variables across all 
sites, but GEMI and the NDVIs from spring (4.03%) and winter (3.93) are 
considered more important than the highest ranked S-2 index IRECI 
(3.44%). Seasonality information from seasonal NDVIs shows again the 
highest values in La Campana 2 (15.16%), while the importance inside the 
vegetation classes is much lower for La Campana 1 (6.43% compared to 
12.25% with all classes). Percentage in Nahuelbuta is higher (11.88%), with 
all seasonal NDVIs showing up in the top 15. 

An exemplary result for class-wise VI is shown for Araucaria in Figure 24. 
The five DEM features together account for 29.28% of the total predictive 
power, which is more than twice as much as within all vegetation classes. All 
seasonal NDVIs appear in the ranking as well as four bands from the infrared 
region. Bands from the infrared region are dominant within the Nahuelbuta 
prediction. Classes 2, 5 and 7 each feature two bands from the infrared re¬ 
gion inside the top 3 most important variables, while only Araucaria feature 
two DEM derivates. In La Campana 1, the DEM is more significant for all 
vegetation classes being for 5 out of 7 classes the most important variable. 
Shrubland alone does not feature a DEM feature in its top 3. In contrast to 
the coastal part, in La Campana 2 DEM is indeed the most important feature 
but for single vegetation classes winter and spring NDVI are more promi¬ 
nent. They appear in every top 3 except for Shrubland. Class-wise analysis 
also confirms that radar features are rated unimportant for every vegetation 
class. Yet, in La Campana 1 GLCM features especially Mean and Variance, 
are under the top 10 most important variables. This is different in La Campa¬ 
na 2, but similar in Nahuelbuta. 


5.6 Time Series Analysis 

Analysis of seasonal NDVI patterns was conducted for the study areas with 
highest vegetation cover (La Campana, Nahuelbuta). Figure 25 plots all 160 
NDVIs after removal of files heavily biased by atmospheric conditions (e.g. 
clouds, fog). However, this does not imply that atmospheric and topographic 
effects such as shadowing and pixel value variation due to slightly different 
illumination conditions could be fully excluded. NDVIs derived from Land- 
sat 8 are in average 0.09 higher than those derived from Landsat 7 or Landsat 
5 (Landsat 8 overall mean: 0.52; Landsat5&7: 0.43). There is no such a clear 
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Figure 25: Complete time series plotted in La Campana (06.08.1997 - 29.08.2017). 

Note the differences in mean NDVI (catchment-wide) between Landsat 8 
and Landsat5&7. Values between data points are interpolated using line¬ 
ar interpolation. Interpolation is disrupted if time gaps between subse¬ 
quent NDVIs is greater than 2 months. 


difference between Landsat 5 and Landsat 7. Due to the failure of Landsat 7 
in mid-2003, several larger gaps exist till the launch of Landsat 8. Between 
12.2011 and 03.2013 no Landsat imagery are available. For the time after the 
launch of Landsat 8, dense time series are available and monthly NDVI has 
just 4 gaps for both La Campana and Nahuelbuta. 

Results of the time series analysis are split into catchment-wide analysis, 
differentiation of NDVI series for different altitudinal belts as well as single 
classes. For the latter a subset of each four vegetation classes is presented. 
BFAST was applied to all kinds of time series but only results for the catch¬ 
ment-wide analysis are depicted as the results of class-wise and height spe¬ 
cific analyses are largely redundant. Some of the most important findings 
from class-wise and height specific analyses are mentioned in the text. The 
reader is referred to the Appendix for BFAST graphs of all time series. For 
the analyses, only data from Landsat 8 was used due to the offset from NDVI 
values of previous Landsat missions and only NDVI >0 were considered to 
filter out shadows, water, ice and snow. 
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5.6.1 Catchment-wide Analysis 

Landsat 8 time series NDVI calculated as mean over the whole catchment are 
shown in Figure 26 and Figure 27. The seasonal variation is visible in both 
study regions (note that seasons on the southern hemisphere are shifted by 6 
months). Flowever, there are differences in the timing of green-up and peak 
vegetation. Peak vegetation is the maximum NDVI after which the slope of 
the NDVI curve changes its sign (Pettorelli et al., 2005). In La Campana the 
date of maximum NDVI is in the summer months, while in Nahuelbuta it is 
early in the year. Green-up phases when defined as the time from minimum 
to maximum NDVI are shorter in Nahuelbuta, while the NDVI stays longer 
on a high level. The green-up in La Campana is longer than in Nahuelbuta 
with a more pronounced peak. Values in the vegetative season of La Campa¬ 
na are around 0.6 and approx. 0.2 higher in Nahuelbuta. 

Analysis of trends with BFAST in R requires regular frequencies, whereby 
the start time of the cycle is unimportant as long as a cycle is defined by the 
same duration (12 months in this case). Therefore, the time series were 
clipped to 48 months between 04.2013 and 03.2017 to minimize gaps. After 
clipping, 4 gaps remained in the data of both study areas. Missing values 
were imputed using spline interpolation. 

BFAST splits the time series data into its seasonal and trend component. 
Figure 28 shows the results of the BFAST analysis for La Campana. The 
fitted seasonal component starts from a local minimum in spring with a 
strong increase in NDVI. The curve flattens out when approaching the date 
of maximum NDVI. Then, a strong decline takes place until the curve flat¬ 
tens out towards the date of minimum NDVI in late spring. The trend com¬ 
ponent identifies a disruption of this cycle in mid-2015 till mid-2016. Be¬ 
forehand, the trend is marked by a slight decline. When the dismption takes 
place, NDVI first falls below its usual level then increases to a level far 
above the normal state. The strength of the decrease marked by BFAST is 
about 0.08, from 0.50 to 0.42, the strength of the increase is about 0.15, from 
0.42 to 0.57. The maximum value in the trend component is 0.63, which is 
slightly lower than the maximum of the original data (0.67). The timing of 
the break within the 95% confidence interval is between April 2015 (25 th 
observation) and April 2016 (37 th observation) and thus the disruption lasts 
almost exactly one year. 
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Figure 26: Zoom-in into the time series from 12.04.2013 - 29.08.2017. The graph 
shows seasonal variation of up to 0.3 in mean NDVI. Highest NDVI val¬ 
ues (peak vegetation, marked with arrows) are between June and Sep¬ 
tember (south winter-spring), while during December to April (austral 
summer-autumn) NDVI is low. Green-up phase preceding the date of 
maximum NDVI is approximated with green shading. 


0.95 



Figure 27: Complete time series for Nahuelbuta from 12.04.2013 - 25.05.2017. The 
graph shows seasonal variation up to 0.4 in mean NDVI. Highest NDVI 
values (peak vegetation, marked with arrows) are between December and 
March (austral summer-autumn), while during June to September (south 
winter and spring) NDVI is low. Green-up phase preceding the date of 
maximum NDVI is approximated with green shading. 
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Figure 28: Analysis of the seasonal patterns and trends with BFAST in the Ocoa 
Catchment (La Campana). Time series from 04.2013 - 03.2017. The ab¬ 
breviations on the y-axis stand for: Yt = is the NDVI (the univariate time 
series under analysis) plotted over time; St = is the fitted seasonal com¬ 
ponent of Yt. St is fixed between breakpoints but can vary across them 
(disturbances); Tt = is the trend component showing disturbances 
(breaks) with dashed lines and the 95% confidence interval for the timing 
of the breaks; et = is the remainder meaning the estimated noise or varia¬ 
tion that remains in the data beyond the seasonal and trend components. 



Time 

Figure 29: Analysis of the seasonal patterns and trends with BFAST in the Nahuel- 
buta catchment. Time series from 04.2013 - 03.2017. No significant 
breaks are detected. The linear equation of the trend component (Tt) is 
given showing a slight increase. 
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The seasonal component in Nahuelbuta inhibits no skewness as found in La 
Campana. The Increase towards the date of maximum NDVI as well as the 
decline afterwards have approximately the same slope. After the start of the 
green-up phase between June to September, the NDVI curve culminates in a 
flat summit. The extreme values of the original data were 0.53 and 0.85. The 
fitted trend component shows a small linear increase, which can be described 
by the formula y = 0.0004x + 0.7169 beginning with a value of -0.72 and 
resulting in -0.74. The trend component as found by BFAST maintains the 
mean of the original data (0.73). Consequently, this is also the case in La 
Campana (0.52), while the seasonal component detected by the algorithm 
always sums up to 0. No disruption is found within the overall trend across 
years. 

5.6.2 Height-specific Analysis 

Based on the TanDEM-X DEM, the study areas were split into equidistant 
altitudinal belts of 200 m. The same spacing was chosen for both La Campa¬ 
na (Ocoa catchment) and Nahuelbuta resulting in 12 different belts for La 
Campana and 5 belts for Nahuelbuta. Areal coverage of single belts is vary¬ 
ing strongly as can be seen within Figure 30 and Figure 31. 

In La Campana (Minimum elevation: 379 m; Max.: 2229 m) the heights 
between 400-600 m feature the highest inter-seasonal variation (Min.: 0.36; 
Max.: 0.76). The mean NDVI across all years in this height zone is 0.57 and 
therefore only slightly lower than the mean in heights between 600-800 m 
(Mean: 0.58), which exhibit the strongest vegetation signal with slightly less 
seasonal variation (Min.: 0.39; Max.: 0.75). From >400 m till 1400-1600 m 
the difference between Minimum and Maximum decreases from 0.51 to 0.24. 
After this point no clear trend is visible, while the lowest variation occurs 
between 2000-2200 m (Diff: 0.21). Standard deviation (SD) shows similar 
behavior with values decreasing steadily from 0.10 (>400 m) to 0.05 (1600- 
1800 m). Note that the seasonal influence is increasing again for heights 
>2200m (Diff: 0.29; SD.: 0.070) being also visually more pronounced, while 
the areal extent is limited to a small area (0.028 km 2 ) close to the summit of 
El Roble (which itself is outside the catchment). The BFAST analyses 
marked significant breaks inside each altitudinal belt except for the zone 
<400 m. 
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Figure 30: NDVI time series by height a.s.l. in the Ocoa catchment (La Campana). 

With time on x-axis and NDVI on y-axis. The small picture shows the 
catchment and the area corresponding to the respective height. The tripar¬ 
tite seasonal component is shown in the time series between 1400-1600 
m. Only a subset of the 12 200 m altitudinal belts is shown. Note the dif¬ 
ferent value range compared to Figure 31. 


Breaks detected in the zones 1400-1600 m and >2200 m correspond well 
with the breaks detected in the catchment-wide analyses, while the seasonal 
components fitted in the zones 800-1000 m, 1000-1200 m and >2200 m are 
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Figure 31: NDVI time series by height a.s.l. in Nahuelbuta. With time on x-axis and 
NDVI on y-axis. The small picture shows the catchment and the area cor¬ 
responding to the respective height. The plot for heights between 800- 
1000 m was omitted because the distribution was almost identical to the 
>800 m plot. 

closest to the seasonal component fitted to the catchment-wide time series. 
The altitudinal belt between 1400-1600 m differs from the others because it 
is the only time series which cannot be described by a single seasonal com¬ 
ponent due to its irregular complexity (Figure 30). 
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The altitudinal gradient in Nahuelbuta is lower than within the Ocoa catch¬ 
ment (Min.: 732 m, Max.: 1534 m). The maximum NDVI in Nahuelbuta is 
0.85 (22.01.2016; 1000-1200 m), the minimum is 0.36 (12.08.2014; 1400- 
1600m). Seasonality is existing in every altitudinal belt, with the highest 
difference being in the highest elevation zone (Diff: 0.46) and the lowest 
between 800-1000 m (0.28). Regarding variation expressed by standard de¬ 
viation, the 1000-1200 m (SD: 0.090) and 1200-1400 m (SD: 0.088) zones 
have higher values than in the belt >1400 m (SD: 0.086). Seasons are more 
pronounced between 1000-1400 m, which covers also most of the area inside 
the catchment (24.32 km 2 ), while apart from two exceptions variation >1400 
m is low. The first outlier is located on the 12.08.2014, the second on 
22.05.2016 after removal of both SD was reduced to 0.057, which is the 
minimum across all altitudinal belts. Highest and lowest elevation zones 
cover only a small percentage of the whole catchment. In the north eastern 
depression, a stream leaves the catchment to ultimately join the Rio Pi- 
coiquen approx, half a kilometer outside the catchment borders (area <800 m 
= 0.40 km 2 ). On the slope of the Nahuel mountain (north-west), elevations 
>1400 cover an area of 0.54 km 2 . The frequency of vegetation classes in this 
area is as follows: Young Growth (1864 m 2 ), Deciduous Forest (1255 m 2 ) 
and Araucaria (779 m 2 ). No breaks are detected within altitudinal belts 
>1000 m and single breaks inside the two time series <1000 m centered on 
April 2016. Seasonal components fitted to the time series are more conform 
than those fitted in different height zones in La Campana. Positive trends are 
found in the height zones 1000-1200 m and >1400. The trend components in 
height zones <1000 m instead show a negative trend before and after the 
break. 

5.6.3 Class-wise Analysis 

Four exemplary vegetation class NDVI time series are depicted in Figure 32 
and Figure 33. In the Ocoa catchment, the highest mean NDVI of all classes 
is 0.80, highest maximum NDVI value is 0.88 (22.05.2016) both within Or¬ 
chards, which are spatially limited to low altitude areas in the catchment's 
northern outlet. Orchards also show little seasonal variation, while the SD is 
higher (0.071) than that of classes 2 and 6 (0.060) due to two outliers in the 
Orchard NDVI time series on 24.12.2013 (0.48) and 06.03.2016 (0.59). After 
removal of these two outliers, SD of class 4 decreases to 0.042. Shrublands 
exhibit the lowest NDVI (0.35; 06.03.2016) and the lowest mean NDVI 
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Figure 32: NDVI time series by vegetation class in the Ocoa catchment (La Campa- 
na). With time on x-axis and NDVI on y-axis. The small picture shows 
the catchment and the area corresponding to the respective class (extract¬ 
ed from the RF classification, visualized based on DEM). Note the edges 
introduced by differences in classification between La Campana 1 and La 
Campana 2 (e.g. Shrubland). 

(0.49). Time series of classes 1, 2 and 3 follow similar patterns of NDVI. 
Simple correlation using the Pearson coefficient shows high values between 
classes 1 and 3 (0.95) as well as 2 and 3 (0.87), while correlation between 
classes 1 and 4 is much lower (0.67). The BFAST analysis detects breaks in 
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Figure 33: NDVI time series by vegetation class in Nahuelbuta. With time on x-axis 
and NDVI on y-axis. The small picture shows the catchment and the area 
corresponding to the respective class (extracted from the RF classifica¬ 
tion, visualized based on DEM). 


each class except for class 4. The timing of the breaks corresponds well with 
the breaks found in the catchment-wide analysis within classes 1, 2, 5 and 7. 
Classes 3 and 5 deviate from the aforementioned due a shift towards 2014 
(1 st break mid-2014, 2 nd break end-2015). The seasonal component of class 3 
is closest to the seasonal component fitted to the catchment-wide analysis. 
Classes 1, 2 and 5 are marked by a longer duration of high NDVIs in the 
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summer months. Mean NDVIs are higher in Nahuelbuta than within vegeta¬ 
tion classes in the Ocoa catchment (Figure 33). Three classes exceed a mean 
NDVI of 0.8 (classes 5, 6 and 12), while the lowest mean values are associ¬ 
ated with the mixture class 14 (0.59) and Grasslands (0.69). SDs range from 
0.04 (classes 5, 6 and 12) to 0.10 (class 2). Likewise, the differences between 
minimum and maximum NDVIs are small for classes 5, 6 as well as 12 
(Min.: 0.75; Max.: 0.91) and greater for class 2 (Min.: 0.55; Max.: 0.88). 
Grasslands and deciduous forests exhibit a trend towards higher NDVIs, with 
a less strong decline in the winter months (2015, 2016) and higher maxima in 
the summer months. These trends are not visible in Araucaria and young 
growth time series. The number of breaks detected in the class specific time 
series with BFAST is higher than within altitude specific time series. Two 
breaks are detected in the time series of classes 2 and 5 and one break is 
detected in the time series of class 12. Disruption in time series of classes 2 
and 5 lasts from early to mid-2014 till mid to late 2015. The single break in 
the time series of class 12 is located at approx, the central position between 
breaks in classes 2 and 5 (early 2015). All other classes (classes 6, 7, 13 and 
14) do not exhibit breaks but a positive trend, which is strongest within clas¬ 
ses 7 and 14. The seasonal components of classes 2 and 13 are closest to the 
seasonal component fitted to the catchment-wide analysis. Classes 5, 6, 12 
and 14 feature a seasonal component with two maxima (with the 1 st being 
higher than the 2 nd ) and a shorter stay on low NDVI values. 


5.7 Summary & Interpretation of Results 

Apart from the focal results (vegetation maps, BFAST graphs), this study 
features a lot of preliminary results of which the most significant are also 
described in this chapter. The initial comparison of topographic corrections 
was the first step that provided insights into the performance of different 
correction algorithms in the La Campana region. Overall, it was shown that 
bands are subject to varying strength of correction. E.g. the deviation from 
the mean was generally stronger in the infrared and less pronounced in the 
Red region of the visible spectrum. Most algorithms also had in common that 
the mean was shifted to lower pixel values leading to an image darkening. 
According to the assumption that lower deviation from and maintenance of 
the mean can be regarded as “better”, Sen2Cor performed worse than any 
one of the other algorithms. Nonetheless, the algorithm did comparably well 
in flattening topography and dealing with different illumination conditions 
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especially inside the mountainous region of La Campana 2. Although some 
overcorrection expressed by bright spots in high topography areas and persis¬ 
tent illumination effects in La Campana 1 existed, Sen2Cor visually outper¬ 
formed simpler algorithms such as the Improved Cosine or Gamma correc¬ 
tion. Regarding the introduction of artifacts on shadowed spots in high to¬ 
pography areas, Sen2Cor showed moderate performances. While some algo¬ 
rithms such as Minslope and C-factor introduced stronger artifacts with 
strange coloring from blue to greenish, others such as Improved Cosine did 
not introduce any artifacts at all. However, in general inside this group of 
algorithms to which also Gamma correction belongs, undercorrection is evi¬ 
dent and the flattening effect is poor. 

Sampling in the difference area showed low accuracy statistics for both at¬ 
mospherically or topographically corrected data. This can be explained by 
the way samples for accuracy assessment were taken. The difference area 
between the two classifications represents mostly pixels, which cannot be 
easily assigned to a specific class i.e. mixed pixels. The error rate at these 
places is certainly higher than elsewhere in the target region. The assessment 
showed that atmospheric correction did overall better as testified through the 
difference area sampling and the classifications based on different prepro¬ 
cessing applied to the data. E.g. an expected increase in classification accura¬ 
cy in the low illumination zone for topographically corrected data did not 
happen. To the contrary solely atmospherically corrected data performed 
better. Albeit, McNemar test demonstrated that the differences are mostly 
insignificant and the ratio of class-wise statistics per preprocessing after 
exclusion of differences <0.2 presented in chapter 5.3.2 even is in favor of 
the topographically preprocessed data. 

The comparison of classifiers showed that RF performed better than SVM, 
although the differences found were only in one case statistically significant. 
None of the two machine learning approaches applied here did better in pre¬ 
dicting each class at a certain study region. Mostly these differences were 
small and only in some cases such as agriculture in La Campana 2 the differ¬ 
ences were evident. In this specific case SVM failed completely in correctly 
predicting agriculture at the validation sites, whereby agriculture was existent 
and predicted correctly at other sites as visual comparison suggests. In com¬ 
parison to SVM, RF has the advantage that VI estimates can be derived, 
which is apart from the slightly higher accuracy statistics another reason why 
RF was used to build the final LULC maps. The OA statistics of the final 
LULC maps are high especially for Pan De Azucar and Santa Gracia. This is 
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due to the relative dominance of bare ground, which covers most of the sur¬ 
face area. The number of samples for training and validation of vegetation 
classes was comparably low, which explains the almost perfect accuracy 
even if misclassification is present. It is therefore advisable to consider the 
class-wise statistics. Vegetation mapping of classes such as shrubs&cacti in 
the desert is subject to the subpixel level. So, pixels assigned to this class in 
both Pan De Aziicar and Santa Gracia are unlikely to exhibit 100% vegeta¬ 
tion cover but rather a certain fraction that allows for a vegetation signal to 
be present. This lead to an overestimation of classes with low fractional veg¬ 
etation cover, in contrast to classes with higher cover such as trees&shrubs 
and irrigation agriculture in the arid study regions. In Santa Gracia shrub 
communities increasingly occur within a certain distance to the coast (100- 
1000 m) at altitudes somewhat over 200 m. This cumulated occurrence of 
low vegetation near the coast can be interpreted as members of the Lomas 
formation. 

A similar dependency of vegetation with location was found in La Campana, 
where exposition determines the growth structure and composition of natural 
vegetation. On south-facing slopes, higher vegetation in form of Mediterra¬ 
nean trees is present, while shrub communities dominate slopes exposed to 
the north. This dependency is conditioned by illumination and its conse¬ 
quences on the local water balance. South-facing slopes are wetter than 
north-facing slopes, which receive more light and loose more water due to 
evapotranspiration. This characteristic distribution of vegetation is present in 
both scenes from La Campana. The overall classification differences between 
the scenes are most obvious in terms of classes 1, 3 and 9. Classes 1 and 9 
substitute for class 3 in La Campana 2, so the RF classifier seems to be more 
likely to classify those. Due to the fluent transition between e.g. bare ground 
and shrubland already small differences in class allocation can lead to major 
differences in the overall classification. 

The analysis of VI gives important insights into the significance of single 
variables from the whole classification process down to single classes. VI of 
a variable allows for valuable assumptions on the data. A high VI inside a 
certain class e.g. suggests that this variable is important to distinguish this 
class from others. Therefore, a high overall VI across all classes indicates 
that classes differ in terms of characteristic values of this variable and that an 
omission would lead to a major reduction in discriminability, which would 
thus reduce classification accuracy. On the contrary a low VI suggests that a 
variable is no characteristic of a certain class, which would allow differentia- 
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tion from other classes. This was the case in Pan De Aziicar where DEM 
features in general had low importance scores, which implies that vegetation 
in this study site does not depend on elevation or topography in general. In 
all other study areas elevation was significantly contributing to classification 
accuracy being the feature with the overall highest MDA. Height was an 
important predictor for both vegetation and non-vegetation classes, e.g. 
snow&ice, shadows as well as deciduous trees in La Campana only occur in 
certain altitudinal belts. It confirms expectations that the DEM was the most 
important feature for each of these classes. Although elevation had the high¬ 
est mean among all vegetation classes in La Campana 2, seasonal NDVI 
features (winter & spring NDVI) did appear more often than DEM features. 
Even in comparison to Nahuelbuta, which features extensive forests of de¬ 
ciduous trees, the relative importance of seasonal NDVIs in La Campana 2 
was larger. As this significance of seasonality features for the classification 
was missing in La Campana 1 it is likely that the higher elevation in La 
Campana 2 is the reason for a stronger influence of seasonal vegetation fluc¬ 
tuations. Radar features had generally low VI, except for the mapping of 
artificial surfaces where GLCM features were considered important. This is 
because in both GLCM Mean and Variance artificial surfaces have high val¬ 
ues, while spectrally similar classes such as agriculture or bare ground had 
low values. Apart from the examples presented above, characteristics of the 
spectral signature can also lead to strong consideration in the classification 
process. In Figure 34, five exemplary vegetation classes from Nahuelbuta 
and their respective mean reflectance values in each S-2 band are plotted. In 
terms of deciduous trees, SWIRI is considered to be of high importance. The 
graph shows that the reflectance values in the SWIR at this place differ from 
those of other classes, so that a differentiation based on this specific Sentinel 
band may indeed be possible. Furthermore, the reflectance in the near infra¬ 
red region differs between single vegetation classes. It is therefore not sur¬ 
prising that in Nahuelbuta as well as La Campana the infrared region espe¬ 
cially the Red Edge 4 band situated at the point of highest reflectance in most 
cases, is a significant predictor for vegetation classes. 

Dynamics and trends in vegetation on catchment scale were assessed using 
NDVI time series and analysis with BFAST. The algorithm detected a break 
in the Ocoa catchment time series lasting from 2015 till mid-2016. Climatic 
forcing is likely to be the reason for this anomalous behavior of the NDVI as 
the Southern-Oscillation Index (SOI) indicated an El Nino event that falls 
into this period (cp. www.bom.gov.au). While the influence of this event 
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Figure 34: Spectral signatures of 5 vegetation classes from Nahuelbuta as seen by 
S-2 with the order of the 10 optical bands above (“RE” = Red Edge). 
Reflectance values calculated as means per band from the training 
samples. 


seems to vanish when entering the temperate climate zone of the Nahuelbuta 
region. Inside the National Park, no break and anomalous behavior of the 
NDVI could be observed. Over time a slight increase in NDVI was observed 
in this catchment. These trends are more pronounced within single classes 
such as deciduous trees and grassland. On the class level, breaks were detect¬ 
ed, but the signal of single classes is attenuated within the overall study area 
mean. The presence and relative dominance of certain classes inside different 
altitudinal belts leads to a strong imprint in the respective height zone time 
series graph. E.g. in heights between 1000-1200 m deciduous trees are domi¬ 
nant in Nahuelbuta, therefore both height and class specific graphs show a 
similar pattern. The same is true for elevations >2200 m in La Campana, 
where again deciduous trees lead to a stronger imprint of seasonality in the 
height specific plot. 
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The main research objectives of this study are the provision of comprehen¬ 
sive vegetation maps combined with the monitoring of vegetation dynamics 
in four study sites along a climatic gradient in Chile. From these focal objec¬ 
tives three different hypotheses are formulated. The discussion chapter is 
divided into three parts each attributed to one of the hypotheses. 


6.1 Hypothesis 1 

Hypothesis 1 states that the vegetation in Chile is distributed heteroge¬ 
neously over the different eco zones along a climatic North-South gradi¬ 
ent, but also along a West-East gradient dependent on the altitudinal 
belts and distance from the coast. 

To answer this question maps of LULC are derived in two arid (Pan De 
Azucar, Santa Gracia), a Mediterranean (La Campana) and a humid- 
temperate (Nahuelbuta) study site. Analysis of these maps supports the pre¬ 
sented hypothesis. The density of vegetation is increasing from the north¬ 
ernmost study site Pan De Azucar where plant coverage is less than 1%, 
towards the south and is culminating in Nahuelbuta, which is by more than 
50% covered with vegetation. The number and sort of vegetation species and 
communities that could be differentiated within the 10 m resolution of the 
Sentinels is varying along with this gradient. While in Pan De Azucar shrubs 
and cacti fonn the major part of species diversity, trees begin to occur in 
Santa Gracia and grow in abundancy towards more favorable climates in the 
south. But there is no steady increase from north to south as found for plant 
cover. The number of vegetation classes in La Campana (7) is higher than the 
number of classes in Nahuelbuta (6). While this alone is no valid factor for 
measuring species diversity of a study site as it largely depends on the ability 
of the sensor to distinguish certain vegetation communities, a number of 
studies support this observation. Various authors beginning from Rundel and 
Weisser (1975) to Hauenstein (2012) and Moreira-Munoz (2011) in his book 
on plant geography of Chile point out on the unique position the La Campana 
region and the National Park have inside the Chilean eco zones. The Medi- 
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Figure 35: Araucaria Araucana stands and outcropping bedrock on a highly elevated 
point inside Nahuelbuta National Park (with kindly permission by © Juli¬ 
ana Klein (2018). 

terranean environment features a great diversity of higher plants ranging 
from shrubs&cacti communities to different forest types thus posing the 
biggest challenge for vegetation mapping. Vegetation in all study sites also 
features a dependency with elevation that is reflected in the importance of the 
elevation variable in RF classifications and the relative abundancy of classes 
depicted in the cumulative plots in chapter 0. Some of the most prominent 
examples for this relationship are the deciduous forest stands of Nothofagus 
obliqua on high elevations in the La Campana National Park (Rundel and 
Weisser, 1975) and the pure Araucaria stands limited to the highest points 
with outcroppings of bedrock (Endlicher and Mackel, 1985) (Figure 35). 

In contrast to this, the relationship with distance from the coast is less obvi¬ 
ous in most of the study sites. The long history of human occupation has 
altered the natural landscape and affected the Chilean vegetation belts 
through localized mining in the north to extended forest substitutions in the 
south (Moreira-Munoz, 2011). The degree of human interference is unequal¬ 
ly distributed across the landscape but focusses on some areas, e.g. major 
settlements in the La Campana region are in close proximity to the coast or at 
least on the valley floor (Valparaiso, Limache, Quillota). All study sites have 
experienced anthropogenic influence, even inside the absolute desert of Pan 
De Azucar roads are existent. It is more difficult to estimate if a certain 
community is more frequent in an area with strong human impact because a 
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major part of the natural landscape has been converted. On the other hand, a 
seemingly existent dependency with distance from the coast may actually be 
a dependency with another variable as it is the case with Araucaria that occur 
within 30-45 km from the coast only because the Cordillera de Nahuelbuta is 
located here. More obvious is the coastal accumulation of low shrubs&cacti 
in the Santa Gracia region, which are interpreted as members of the Lomas 
formation. 

The LULC maps are produced by a RF classification framework that outper¬ 
formed SVMs both in terms of OA and in most cases also class-related accu¬ 
racies. RFs were successfully applied in various LULC studies before and 
most often compared to simple classifiers such as Decision Trees or other 
ensemble learners such as bagging and boosting (Gislason et al., 2006; 
Waske and Braun, 2009; Rodriguez-Galiano et al., 2012). Comparative stud¬ 
ies with other approaches from the field of machine learning are rather rare. 
One of the few studies was conducted by Pelletier et al. (2016), who used 
SPOT-4 and Landsat-8 high-resolution time series to map LULC in southern 
France. Their comparative classification approach indicated a 6% higher OA 
for RF (83.3%) than for SVM (77.1%), whereby the RF classifier also 
showed a couple of other advantages including lower computational cost and 
easy parametrization. Indeed, the present study observed a up to two times 
faster performance of RF for both parametrization and prediction. In addi¬ 
tion, RF has the advantage that variable importance can be calculated. In 
contrast to the black-box approach of SVMs, this allows at least for some 
degree of transparency, feature selection and the estimation of significance 
for the whole classification as well as single classes. Overall the quality of 
the produced LULC maps is satisfying and the results confirm that S-2 has 
high potential for land cover mapping (Immitzer et al., 2016). The achieved 
accuracies are on the same level than those reported by Borras et al. (2017) 
in their S-2 land cover mapping in Argentina and Spain with a similar num¬ 
ber and composition of classes. In addition, all maps exceed the accuracies 
reported by Zhao et al. (2016) in their study across all levels of detail. While 
their number of classes distinguished within level 1 is 12, this study has 19 
classes although not all are present inside a single study site. The least clas¬ 
ses could be distinguished within the desert environment of Pan De Azucar. 
LULC mapping in arid environments is generally hampered by the low vege¬ 
tation cover and sand covering anthropogenic features such as roads. Inside 
the Atacama Desert, the vegetation if existent mainly consists of withered 
shrubs and small cacti. Even within higher resolution imagery as available 
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Figure 36: Difficult to detect within 10 m resolution: Sparse shrubs close to Pan De 
Azucar National Park (with kindly permission by © Juliana Klein (2018). 

through GoogleEarth, vegetation is hard to identify and easily confused with 
small rocks. 

Figure 36 shows shmbs close to the Pan De Azucar National Park. The di¬ 
ameter of an individual is approx. 40 cm. In addition, most plants especially 
those further away from the coast are dryer and lack a green coloring that 
could help in their recognition. The reported accuracy statistics in Pan De 
Azucar and Santa Gracia are overestimations. The misclassification especial¬ 
ly in terms of shrubs&cacti is likely to be significantly higher. Separability 
measures also indicate potential issues in distinguishing shrubs and bare 
ground through a value of 0.37, which can be interpreted as “poor separabil¬ 
ity”. Only pixels that could be clearly assigned to a vegetation class were 
sampled, which lead to much less samples for training and validation com¬ 
pared to other classes especially bare ground. Indeed, these pixels with a 
clear vegetation signal were mostly assigned to the correct class, but pixels 
with a higher contribution of bare ground and rocks of colors similar to 
plants were often falsely attributed. This leads to an overestimation of the 
area in which low vegetation coverage classes such as shrubs&cacti occur. 
The same is reported by Zhao et al. (2016) where the OA in the Atacama 
Desert region was the lowest among all ecoregions due to the confusion of 
withered shrubs with barren lands. Another factor that must be taken into 
consideration is the approach chosen for sampling in this study. Despite the 
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Figure 37: Do two palm trees in a 10 m 2 plot make a “Mediterranean forest” pixel? 

Example from the Ocoa catchment (with kindly permission by © Juliana 

Klein (2018). 

sampling was done carefully to systematically cover the whole image and the 
whole spectral range as far as possible while sampling no pixel more than 
once, the non-random approach introduces some autocorrelation and nega¬ 
tively affects the validity of the accuracy measures. However, as discussed in 
chapter 4.1.3 the non-random approach that samples only at sites, which can 
be clearly assigned to a certain class is a compromise to ensure high classifi¬ 
cation accuracy without field validation. This non-random approach was also 
used by Vanonckelen et al. (2013) for training data generation but in combi¬ 
nation with field data. 

A mixed pixel issue between vegetation and bare ground is not indicated by 
the separability measures in the other study sites. In all sites and within all 
class combinations the separability does not fall below 0.5. Albeit, from 
visual analysis and according to the class-wise statistics it is obvious that 
there is confusion between the shrubland and bare ground classes. This is 
again due to the fact that only shrubland with an evident vegetation signal 
was sampled so that the misclassification occurs at pixels of a higher share of 
bare ground. Already small differences in class allocation can lead to major 
changes in the mapping as can be seen in the different scenes of La Campa- 
na. The edges between both classification maps show that in La Campana 1 
the classifier was more likely to assign pixels to shrubland, while in La 
Campana 2 Mediterranean forest and bare ground occur. The transition be¬ 
tween Mediterranean forest and shrubland is fluent as shrubs can play a dom- 
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inant role in both classes especially at north-facing slopes where vegetation 
is dryer and of lower growth. Even specialists from the field of ecology will 
experience difficulties when asked to draw distinct lines between Mediterra¬ 
nean forests and shrub land. An example is shown in Figure 37 from a north¬ 
facing slope in the Ocoa catchment with the Cerro La Campana in the back¬ 
ground. The area is sparsely covered by Palm trees (Jubaea Chilensis ) and 
most of the surface on the slope is covered by shrubs and bare ground. With¬ 
in 10 m 2 plots two palms may cover a substantial part. So, if only a single 
plot is considered without its surroundings it may be assigned to a “Mediter¬ 
ranean forest” class, although this is highly arguable. 

Thus it is not surprising that such a differentiation is difficult within 10-20 m 
as provided by the Sentinels. Stratoulias et al. (2015) noted that the spatial 
resolution of S-2 (in particular the bands with 20 m resolution) could be too 
low to map highly fragmented forests and landscapes. However, this is a 
matter of the scale and resolution intended. The objective of this study does 
not require a cm-scale separation of forest and shrubland communities but 
rather aims for a large-area mapping without delineation of single individu¬ 
als. Within these requirements the achieved accuracies both regarding overall 
and class-wise statistics in all study sites can be considered good with few 
exceptions in classes that are beyond the focus of this study (agriculture and 
shadows). 


6.2 Hypothesis 2 

Hypothesis 2 states that the ENSO phenomenon has a significant impact 
on the vegetation, but the influence is uneven depending on latitude, 
species and altitude. 

In the short time scale that Landsat 8 data are available it is impossible to 
verify the predicted shorter and delayed greening-up season as a consequence 
of climate change in Chiles semi-arid regions (La Maza et al., 2009). Never¬ 
theless, a disturbance of the seasonal trend was detected in the analyzed 
Ocoa catchment (La Campana National Park) that falls into the time the 
2015/16 El Nino affected the climate in Chile. El Nino events in north- 
central Chile (30°-35° S) are characterized by enhanced precipitation during 
the austral winter months (June-August), while in summer the influence is 
considerably weaker (Montecinos and Aceituno, 2003). The first break was 
detected in April 2015, which is quite close to the period in which El Nino 
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influence is at the maximum. To verify this relationship a correlation with 
climate variables (temperature and precipitation) would be necessary. Such 
an analysis could also give an answer to the question why the NDVI curve is 
first falling below before rising above its usual level. As a common El Nino 
event in Chile in the summer months of this region of Chile is characterized 
by high amounts of precipitation that exceed the usual level of humidity in 
the dry regions of the country by several orders of magnitude, the initial 
decrease in NDVI indicating a decline in plant productivity remains unclear. 
Eventually, there is some connection with the overall climatic situation in 
Chile that governed the rainfall household in the years before this event since 
the 2010 La Nina. From 2011-2015, Chile was exposed to a major drought, 
the so called central Chile megadrought (MD). During the MD annual pre¬ 
cipitation deficits ranging from 55-75% from 30° to 38° S were recorded 
causing a major plant productivity deficit in semiarid Chile. Albeit, the ef¬ 
fects on plants were irregular and natural vegetation, e.g. shrublands in the 
semiarid sector showed a strong decline in productivity, while irrigated 
croplands and fast growing plantations of exotic species (Eucalyptus, Pine) 
were less affected (Garreaud et al., 2017). 

In their study on habitat representation of Nothofagus genera after applying a 
dispersal constrained climate change scenario, Alarcon and Cavieres (2015) 
(AlarconandCavieres,E emphasized that species react on changing climate 
each on their own specific way. This implicates that a trend found for a 
whole biome, may not relate to the way especially less dominant elements 
such as understory plants will react. To account for the expected differences 
in behavior of single species and ecosystems, the study sites were split based 
on the classifications to study the pattern of NDVI for each class. The breaks 
found in La Campana were detected in many (classes 1, 2, 5 and 7) but not 
all classes (3, 4 and 5) and height zones. This indicates that the factor(s) that 
lead to this disruption of the NDVI pattern is not affecting each ecosystem 
and species in a similar way. In case of an El Nino event being the driving 
factor this could imply that the climatic consequences accompanying such an 
event lead to an ultimately stronger growth as indicated by high NDVI val¬ 
ues. The MD fits well into this interpretation because the NDVI is increasing 
again after it decreased in the years before 2015 when the drought was most 
intense. 

While natural vegetation classes exhibit a strong relationship to the climate, 
orchards do not show any response and irregularity in the NDVI pattern 
across the years 2015-2016. This is in accordance with the findings of Gar- 
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reaud et al. (2017), who explained this with the irrigation that makes or¬ 
chards independent from natural rainfall. In contrast to their observations. 
Eucalyptus and Pine do behave differently during the drought and after the 
onset of the El Nino. This is because no plantations are existent in the catch¬ 
ment. The mapped occurrences often represent small stands of few individu¬ 
als that found their way into the park without human intervention. These 
individuals are more susceptible to climate variations as they are not growing 
under the favorable conditions of huge plantations. An immediate conse¬ 
quence of the missing response of orchards to climate signals is that also no 
break could be detected in the altitudinal belt <400 m because orchards are 
the prevailing vegetation community. Overall, the results in both catchments 
suggest that species is more important than altitude regarding the disturbing 
effects. E.g. the altitudinal belt >2200 m in the Ocoa catchment resembles 
that of the deciduous forest class showing the same timing of the breaks 
much because it is governed by the deciduous Nothofagus obliqua forest 
stands. To verify this assumption and exclude that the NDVI curve of species 
varies with altitude an additional combined altitude & species-specific analy¬ 
sis would be necessary. 

The NDVIs in the humid temperate Nahuelbuta catchment are generally 
higher than those in the Mediterranean Ocoa catchment. Considering NDVI 
values >0.8 one has to take into account that the NDVI tends to saturate with 
high vegetation cover and biomass (Gill et al., 2009). This may lead to an 
underestimation of the actual amount of vegetation present during the sum¬ 
mer-autumn season of maximum NDVI. In contrast to the Ocoa catchment, 
no breaks were detected within the Nahuelbuta catchment. If the time series 
are split into class/altitude specific time series, the BFAST algorithm indicat¬ 
ed breaks in some cases (especially classes 2 and 5) but none of the found 
disturbances equals those found in the Ocoa catchment. Furthermore, the 
start of the disturbances (in case two breaks are found ) is already in 2014 and 
thus before the onset of the 2015/16 El Nino. Instead, a gradual increase is 
present in those classes without breaks (classes 7, 13 and 14) and after the 
breaks in the other classes. It seems plausible that this represents a recovery 
of the vegetation from the impacts of the MD. After long years of reduced 
rainfall and increased evopotranspiration (Garreaud et al., 2017), vegetation 
growth and primary productivity are increasing again. For the long term, 
climate predictions suggest a decrease in annual rainfall in the Araucaria 
region of as much as 19-43% depending on the scenario for representing CO 2 
emissions (Orrego et al., 2016). The vegetation response in the region during 
the MD thus has important implications for the future in which such events 
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are likely to happen more frequently (Alarcon and Cavieres, 2015; Garreaud 
et al., 2017). 


6.3 Hypothesis 3 

Hypothesis 3 states that DEM features are valuable auxiliaries in vegeta¬ 
tion mapping. The high-resolution TanDEM-X DEM can compete with 
established global DEMs. 

As pointed out in the discussion of hypothesis 1, the VI measures derived 
from RF support the first part of this hypothesis. Except from Pan De Azucar 
where vegetation is nearly non-existent, the DEM is always one of the most 
important and in case of the VI for vegetation classes only, the most im¬ 
portant feature among all study sites. This further implies that elevation is a 
major control on vegetation inside almost all eco regions of Chile irrespec¬ 
tive of the fundamental differences between those. The relationship with 
elevation of single LULC classes is primarily due to topography influencing 
the patterns of temperature and moisture, which are pivotal location factors 
for plant sustainability and growth (Balzter et al., 2015; Barrett et al., 2016). 
In the Santa Gracia study region, most of the vegetation concentrates on 
lower elevations. Close to the coast humidity through fog favors the presence 
of the extensive shrub and cacti community of the Lomas formation (Rundel 
and Dillon, 1998). The moisturizing effect diminishes further inlands as the 
fog is also hold back by topographic barriers. Higher vegetation and trees of 
the Prosopis genus further inlands exist almost exclusively in topographic 
depressions and on the bottom of valleys. Water from the surrounding 
hillsides and (non-perennial) rivers gathers at these sites thus increasing the 
moisture level and causing more favorable conditions for plant growth. Con¬ 
versely, deciduous Nothofagus obliqua in the La Campana National Park and 
pure Araucaria stands in the Nahuelbuta National Park are both limited to 
high elevations because the conditions (e.g. less competition) are more fa¬ 
vorable for their survival (Rundel and Weisser, 1975; Endlicher and Mackel, 
1985). 

Besides these natural influences of topography on elevation, human land use 
has further contributed to a topography-dependent vegetation pattern in 
Chile. Agricultural irrigation is present in northern semi-arid Chile on the 
border to the Atacama Desert. Especially vineyards have significantly grown 
in spatial extent during the last decade (Young et al., 2010). In contrast to 
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plantations further south in the Mediterranean, vineyards in the main Elqui 
valley and its surroundings focus on the valley floor while the expansion up 
the mountainsides is driven by the increasingly limited space on the valley 
floor. The difference between both study sites can be attributed to insolation 
being no limiting factor in the Santa Gracia region, while the better accessi¬ 
bility, reduced evopotranspiration and better capabilities for irrigation make 
the vineyards on the valley floor more favorable for land owners and compa¬ 
nies. In the La Campana region, deforestation and radiation of anthropogenic 
land use have led to a decrease of natural biomes including Matorral and 
natural forests over time. Evergreen sclerophyllous forest have remained in 
drainage corridors and on steeper slopes with southern aspects. In vicinity to 
rivers, agricultural expansion favored by the availability of water have led to 
the replacement of sclerophyllous forest with Matorral vegetation and elimi¬ 
nation of almost all natural vegetation on the riverbanks (Schulz et al., 2010; 
Vergara et al., 2013). However, in their land cover change study Schulz et al. 
(2011) found that regeneration of shmbland is more likely to occur on flatter 
slopes and closer to rivers. 

Already back in the 70s, Parsons (1976) showed that the land use pressure on 
natural vegetation has altered the characteristics of the Chilean Matorral. He 
compared different species and their adaptation to environmental conditions 
with those of the comparable Californian Chaparral vegetation. Different 
solutions to the same problems (heat, limited water availability) existed with¬ 
in the two study sites. While in the Californian Chaparral deciduous species 
were more common, evergreen species dominated in the more open Chilean 
Matorral. He attributed both the less dense population with species as well as 
the differences in adaptation strategies to a higher degree of disturbance in 
the Chilean Matorral region. These findings become even more important as 
they do have significance for the protected areas, too. Even after being de¬ 
clared a National Park in 1966 illegal charcoal burning and grazing of cattle 
and sheep have persisted in the La Campana park (Rundel and Weisser, 
1975; Hauenstein, 2012). Such disturbances have also been reported from the 
Nahuelbuta National Park (Endlicher and Mackel, 1985). 

From the previous paragraph it has become clear that elevation is an im¬ 
portant factor in vegetation mapping due to interrelationships with tempera¬ 
ture, moisture and human intervention. But also other topographic factors 
derived from the DEM are important indicators for vegetation occurrence. In 
the La Campana region it is most obvious that south-facing slopes differ 
from those facing the equator in terms of the natural vegetation composition. 



6.3 Hypothesis 3 


93 


The wetter southward slopes exhibit sclerophyllous forests, while ridgetops 
and equatorial-facing slopes are dominated by xeric matorral, with spiny 
shrubs and cacti as dominant vegetation (Vergara et al., 2013). Yet, even in 
classes obviously dependent on aspect and sun irradiance in La Campana 1 
(natural forest and shrubland), aspect is ranked with low importance. Instead, 
slope and topographic openness are important topographic variables apart 
from the DEM. Topographic openness describes the dominance or enclosure 
of a pixel, which is related to how far a pixel can be viewed from any posi¬ 
tion (Yokoyama et al., 2002). The openness is positive if the pixel is located 
on a topographically dominant position, which is the case at hillsides and 
ridges. Although VI does not give a hint on how the variable is important for 
a certain class, e.g. if high values or low values are more descriptive, it is 
most likely that positive topographic openness and moderate to high slopes 
are explaining the distribution of native sclerophyllous forests and shrubland. 

Similar to the La Campana region, in the Nahuelbuta study site native forest 
is primarily lost at lower altitudes where accessibility is better and agricul¬ 
ture is more profitable (Altamirano et al., 2013; Miranda et al., 2015). As the 
topography in this region is lower with less steep slopes, the accessibility and 
the vulnerability towards native forest loss is generally greater (Zamorano- 
Elgueta et al., 2012). Native forest therefore remains mainly inside the Cor¬ 
dillera on steep slopes which are less appropriate for agricultural use and 
plantations that dominate on gentle slopes (Endlicher and Mackel, 1985). In 
comparison to the other study sites elevation is less important for vegetation 
classes, although being still the most important variable in case of VI with 
vegetation classes only. This is a consequence of the lower topography in the 
study region, which makes vegetation discrimination based on altitude less 
feasible, except for Araucaria, which are limited to the highest and most 
exposed summits inside the Nahuelbuta National Park. VI is consequently by 
far the most important variable for prediction of this class. 

To assess the performance of the TanDEM-X DEM in topographic correction 
in comparison to other established global DEMs, descriptive statistics were 
derived (Riano et al., 2003) and the difference area between classifications 
based on once a solely atmospherically corrected dataset and once a topo¬ 
graphically corrected dataset were sampled randomly (Vanonckelen et al., 
2013). Compared to Riano et al. (2003) the deviation from the image mean 
and SD were lower in the present study, which suggests an overall good 
performance of the tested correction algorithms. Although visual comparison 
shows artifacts especially in high topography areas introduced by the majori- 
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ty of algorithms including Sen2Cor. Maintenance of means and SD after 
correction was poor for Sen2Cor, while visual correction shows a high po¬ 
tential in flattening moderate terrain and removing topographic effects. A 
number of previous studies also suggest a better performance of Sen2Cor in 
comparison to simple algorithms like Minnaert and C-NDVI (Louis et al., 
2016; Gascon et al., 2017). Thus, the descriptive statistics alone seem to be 
no sufficient measure of the effectiveness of a topographic correction. 

Numerous studies have shown the value of topographic correction including 
Vanonckelen et al. (2013), who found that the impact of the topographic 
correction on classification accuracy was even higher than that of the atmos¬ 
pheric correction. This is of course dependent on the characteristics of the 
study are and satellite images used but shows nonetheless that is worthwhile 
to invest in both atmospheric and topographic correction. Especially in case 
of multi-temporal data as it was the case in the present study. Overall accura¬ 
cies calculated based on the difference area samples showed a higher accura¬ 
cy of the classification based on solely atmospherically corrected data in high 
and low illumination zones, while the topographically corrected data did 
better in the moderately illuminated zone. These results are surprising and 
contradict the expectations. However, the found discrepancies are mostly 
small. If only differences between both classifications >0.2 are considered 
the overall ratio of class-wise accuracies is in support of topographic correc¬ 
tion. All in all, this study could not find a clear benefit but also no general 
drawback from topographic correction. No preprocessing treatment generally 
improved the classification of vegetation classes throughout all illumination 
zones. Of course, this needs further testing on other study sites in order to 
investigate if topographic correction elsewhere increases vegetation classifi¬ 
cation accuracy. Furthermore, the sampling is critical and ground truth would 
ensure the validity of the samples and the results. Finally, other more sophis¬ 
ticated proprietary software such as ATOCR for IMAGINE (GEOSYSTEMS 
GmbH, 2017) could be applied and the results compared to those obtained 
from Sen2Cor. 

Between classifications based on data processed with either SRTM or Tan- 
DEM-X DEM no statistically significant differences could be found. The 
results are slightly in support of the SRTM DEM, but none did show an ob¬ 
viously better performance in predicting vegetation classes in general. The 
TanDEM-X DEM has to deal with issues on water areas (erroneous heights 
and artifacts) and shadowed areas (radar shadowing). In the difference be¬ 
tween Lidar DEM and TanDEM-X DEM presented in chapter 5.3.2, it is 
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visible that the offset is higher on east, south-east facing slopes. This is due 
to the geometry of the acquisition with the radar source coming from the east 
on the descending pass with a right-looking system (Krieger et al., 2007). 
During the SRTM data collection, every point on the earth was recorded at 
least once in ascending and descending pass (Farr et al., 2007), while during 
the TanDEM-X mission phase only high mountainous terrain was mapped 
with both ascending and descending passes (Krieger et al., 2007). The error 
in areas pointing away from the radar beam thus hampers the overall agree¬ 
ment of the TanDEM-X DEM with real-world elevation a.s.l. Additional 
testing would be necessary to investigate if this bias in the TanDEM-X DEM 
hampers the classification within different topography and illumination zones 
since the vertical accuracy and general validity of the TanDEM-X DEM has 
been assured by various studies (Rizzoli et al., 2012; Wessel et al., 2014; 
Balzter et al., 2015). 
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The focal objective of this study was to provide an assessment of state and 
dynamics of Chilean vegetation inside 4 sites along a latitudinal gradient. 
This was accomplished by producing extensive LULC maps using machine 
learning algorithms (RF & SVM) and calculating NDVI time series for a 4- 
year period between 2013 and 2017. The BFAST algorithm (Verbesselt et 
al., 2010) was applied to statistically analyze the time series catchment-wide, 
vegetation category specific and inside altitude belts to investigate the impact 
of climatic forcing. Furthermore, the performance of the recently released 12 
m resolution TanDEM-X DEM in topographic correction and as an input in 
the classification variable framework was tested and compared to a Lidar 
DEM and the established SRTM DEM. 

LULC was predicted with high accuracy (OA: 86-99%). RF did better than 
SVM (2-4%), but differences exist from class to class. The merging of both 
classifiers may lead to more accurate results for some categories than a single 
classifier. Sentinel 10 m resolution data are a major step forward in the open 
availability of high-resolution satellite imagery and the Landsat archive. 
While highly fragmented and marginal plant communities in desert regions 
are still at the limit of being correctly represented. Vegetation in the study 
sites varies not only with latitude but also with distance from the coast and 
topographic factors such as elevation, aspect and slope. Elevation was identi¬ 
fied as the most important variable for predicting vegetation and LULC in 
general throughout all study sites except from Pan De Azucar in the Atacama 
Desert. Against expectations, there is no improvement when using the Tan- 
DEM-X DEM instead of the SRTM DEM in the classification process or 
even if topographic correction is completely omitted. This has been attribut¬ 
ed to effects of shadowing, direction of the radar source and errors above 
water still being present in the DEM. The Sen2Cor algorithm showed deficits 
within strong topography hampering the quality of the corrected products. 

The time series analysis suggests the impact of two major climatic phenome¬ 
na on the vegetation in two Chilean National Parks. A greening-up trend 
from 2013 to 2017 was found in catchments inside both the Mediterranean 
La Campana (Ocoa catchment) as well as the humid-temperate Nahuelbuta 
National Park. It is likely that the increasing vegetation signal detected by the 
NDVI is due to the recovery from a drought occuring in central Chile be- 
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tween 2011 and 2015. A disruption of the NDVI pattern in the Ocoa catch¬ 
ment coincides with the occurrence of the 2015/16 El Nino event. The im¬ 
print of this event is not present in each altitudinal belt or vegetation catego¬ 
ry. This, together with the absence of a disruption in Nahuelbuta, supports 
the conception that the effect of climatic variations on plants in Chile is de¬ 
pending on variables such as latitude, altitude and species. 

Further research should therefore emphasize studies on the small-scale eco¬ 
system and species level rather than observing whole biomes. Time series 
analyses are a valuable source of information about ongoing and past chang¬ 
es as well as trends that may affect future conditions of a species or ecosys¬ 
tem under consideration. In addition, further enhancement is possible by 
correlation with climatic variables such as temperature and rainfall. Although 
the IT-technologies and framework conditions (e.g. Cloud Computing; Open 
Source data) are increasingly made available, there are still less studies, 
which took profit for the purpose of sophisticated time series evaluation. 
Especially in Chile, there is a research gap in terms of satellite time series 
studies with no comprehensive research yet available as it is the case for 
some parts of the Amazonian rainforest. The trends of the last 30 to 40 years, 
in which the archive of satellite imagery has expanded drastically largely 
remain untapped. In view of the rapidly increasing availability of Sentinel 
optical and radar data, there is a huge potential for future research. 
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